diff --git a/CER-calculate-AR.do b/CER-calculate-AR.do index e333468e957209aaf7c362bf5b9b1250f210280e..942c44ac7e093fbe24d849e83a6bd93b08f716d3 100644 --- a/CER-calculate-AR.do +++ b/CER-calculate-AR.do @@ -34,12 +34,16 @@ GNU General Public License for more details. global where "/Users/ben/Documents/Work" global proot "$where/Projects/ESRC-Transformative-Census2022" -global rpath "$proot/results/CER-Irish-SM-Trial" + * original files global odfiles "$where/Data/Social Science Datatsets/CER Smart Metering Project" + * processed files -global pdfiles "$proot/data/cer" +global pdfiles "$proot/data/cer/CER_OctHH_data" + * results path +global rpath "$proot/results/CER-Irish-SM-Trial" + global version "v1" @@ -60,93 +64,90 @@ timer on 1 * see https://github.com/dataknut/Census2022/blob/master/CER-data-processing.do * load 1/2 hour data (don't need survey for this) -use "$pdfiles/CER_OctHH_data/CER_Oct2009HH_30min_survey.dta", clear +use "$pdfiles/CER_Oct2009HH_30min_survey.dta", clear local midweek "0 1" // 0 = weekends -local midweekt1 "Mid-week" +local midweekt1 "Midweek" local midweekt0 "Weekends" +* we have 36 values 'per day' (not 48) +* so lag 35 = this time tomorrow +* for mid-week, this time next week = lag 105 +* for weekend, this time next week = lag 70 +* so set max lags = 150 to cover both (will also mean 2 weekends) +local max_lag = 150 +local midweektl0 xline(35 70, lstyle(refline )) +* the labels don't seem to show up? +local midweektll0 text(35 0.5 "This time tomorrow", place(w)) text(70 0.5 "This time next weekend", place(e)) // weekend -local midweektl0 xline(36 72, lstyle(refline )) -local midweektll0 text(36 0 "This time tomorrow") text(72 0 "This time next weekend") // weekend - -local midweektl1 xline(36 72 108, lstyle(refline )) -local midweektll1 text(36 0 "This time tomorrow") text(72 0 "This time the day after tomorrow") text(108 0 "This time next week") // midweek - -* this creates a big list of all IDs so we can loop over it - takes a while -* could just start from min & loop to max - but would then be testing for non-existent households -qui: levelsof ID, local(ids) +local midweektl1 xline(35 70 105, lstyle(refline )) +* the labels don't seem to show up? +local midweektll1 text(35 0 "This time tomorrow", place(w)) text(70 0 "This time the day after tomorrow", place(w)) text(105 0 "This time next week", place(w)) // midweek local econt1 "In work" local econt2 "Unemployed" local econt3 "Retired" local econt4 "Caring for relative or family" +* filter for +gen s_hour = hh(s_halfhour) + +* this creates a big list of all IDs so we can loop over it - takes a while +* could just start from min & loop to max - but would then be testing for non-existent households +* qui: levelsof ID, local(ids) + +* for testing +local ids "1002 1003 1004" + foreach m of local midweek { di "****************" - di "* Calculating ARs for midweek = `m'" + di "* Calculating ARs for midweek = `m' (`midweekt`m'') " foreach id of local ids { di "* -> ID = `id'" preserve - keep if ID == `id' & midweek == `m' & ds_halfhour >= 12 & ds_halfhour <= 48 // 06:00 - 24:00 -> 36 hour day - * qui: su ds_halfhour - * check - * di "Halfhours: `r(min)' - `r(max)'" - * skip over if fails due to lack of hubid - capture noisily { - * use s_datetime so ac knows what is missing - qui: tsset s_datetime - * do this here so local is set to the unique valuye for this hh - qui: levelsof ba_nchildren, local(nch) - qui: levelsof ba_nadults, local(ba_nadults) - qui: levelsof ba_empl, local(econ) - di "ID: `id' (`midweekt`m'', N = `ba_nadults' adults, `nch' children, respondent: `econt`econ'')" - * we have 36 values 'per day' (not 48) - * so lag 36 = this time tomorrow - * for mid-week, this time next week = lag 108 - * for weekend, this time next week = lag 72 - * so set max lags = 144 to cover both (will also mean 2 weekends) - * this draws graph - slow - qui: ac kwh, gen(archr) lags(150) name(ac_hubid`id'_`midweekt`m'') /// - `midweektl`m''/// draw lines - `midweektll`m'' /// draw line labels - note("ID: `id' (`midweekt`m'', N = `ba_nadults' adults, `nch' children, respondent: `econt`econ'')") - graph export "$rpath/graphs/kwh_archr_hubid`id'_`midweekt`m''_$version.png", replace - - di "* save out arch results for household ID = `id'" - qui: su t - qui: egen lag = fill(1(1)`r(max)') - qui: keep lag archr ID - qui: save "$rpath/tmp/archr-`id'_`midweekt`m''_$version.dta", replace - di "* Done" - } + * take out 00:00 - 06:00 as repeat cycle of sleep will swamp other cycles + keep if ID == `id' & midweek == `m' & s_hour > 6 + * need to fool ac into thinking these observations are continuous + qui: su ID + qui: egen lag = fill(1(1)`r(N)') + tsset lag + * do this here so local is set to the unique valuye for this hh + qui: levelsof ba_nchildren, local(nch) + qui: levelsof ba_nadults, local(ba_nadults) + qui: levelsof ba_empl, local(econ) + qui: su s_hour + di "Halfhours: `r(min)' - `r(max)' for ID: `id' (`midweekt`m'', N = `ba_nadults' adults, `nch' children, respondent: `econt`econ'')" + * this draws graph - slow + qui: ac kwh, gen(archr) lags(`max_lag') name(ac_hubid`id'_`midweekt`m'') /// + `midweektl`m'' /// draw lines + `midweektll`m'' /// draw line labels + xlabel(0(35)140) /// + note("Halfhours: `r(min)' - `r(max)' for ID: `id' (`midweekt`m'', `ba_nadults' adult(s), `nch' children, respondent: `econt`econ'')") + graph export "$rpath/graphs/kwh_archr_hubid`id'_`midweekt`m''_$version.png", replace + + di "* save out arch results for household ID = `id'" + qui: keep lag archr ID + qui: save "$rpath/tmp/archr-`id'_`midweekt`m''_$version.dta", replace + di "* Done" restore } } -* reset the id loop values -qui: su ID -local id_min = `r(min)' -local id_max = `r(max)' -local id = `r(min)' - clear * pool all the results foreach m of local midweek { di "****************" di "* Loading `midweekt`m'' files" - while `id' < `id_max' { - capture noisily { - append using "$rpath/tmp/archr-`id'_`midweekt`m''_$version.dta" - erase "$rpath/tmp/archr-`id'_`midweekt`m''_$version.dta" - } - * distribution of average values for first 70? - graph box archr if lag < 144, over(lag) - graph export "`rpath'/graphs/box_mean_archr_all_hubids_`midweekt`m''_`version'.png", replace - save "`rpath'/archr-all_hubids_`midweekt`m''_`version'.dta", replace + foreach id of local ids { + append using "$rpath/tmp/archr-`id'_`midweekt`m''_$version.dta" + erase "$rpath/tmp/archr-`id'_`midweekt`m''_$version.dta" } + * distribution of average values for first 70? + graph box archr if lag <= 105, over(lag) + graph export "$rpath/graphs/box_mean_archr_all_hubids_`midweekt`m''_$version.png", replace + save "$rpath/archr-all_hubids_`midweekt`m''_$version.dta", replace }