* Fruitfly data -- revised April 3 , to include code for weekly follow-up segments

clear

input id partners  active  lngevity  thorax  sleep

01 8 0 35 0.64 22
02 8 0 37 0.68  9
03 8 0 49 0.68 49
04 8 0 46 0.72  1
05 8 0 63 0.72 23
06 8 0 39 0.76 83
07 8 0 46 0.76 23
08 8 0 56 0.76 15
09 8 0 63 0.76  9
10 8 0 65 0.76 81
11 8 0 56 0.80 12
12 8 0 65 0.80 15
13 8 0 70 0.80 37
14 8 0 63 0.84 24
15 8 0 65 0.84 26
16 8 0 70 0.84 17
17 8 0 77 0.84 14
18 8 0 81 0.84 14
19 8 0 86 0.84  6
20 8 0 70 0.88 25
21 8 0 70 0.88 18
22 8 0 77 0.92 26
23 8 0 77 0.92 24
24 8 0 81 0.92 29
25 8 0 77 0.94 27
01 0 9 40 0.64 18
02 0 9 37 0.70  6
03 0 9 44 0.72 19
04 0 9 47 0.72  7
05 0 9 47 0.72 16
06 0 9 47 0.76 13
07 0 9 68 0.78 35
08 0 9 47 0.80  2
09 0 9 54 0.84 35
10 0 9 61 0.84  6
11 0 9 71 0.84 15
12 0 9 75 0.84 14
13 0 9 89 0.84 18
14 0 9 58 0.88 50
15 0 9 59 0.88 25
16 0 9 62 0.88 10
17 0 9 79 0.88 33
18 0 9 96 0.88 43
19 0 9 58 0.92 35
20 0 9 62 0.92 17
21 0 9 70 0.92 27
22 0 9 72 0.92 22
23 0 9 75 0.92 16
24 0 9 96 0.92 20
25 0 9 75 0.94 37
01 1 0 46 0.64 23
02 1 0 42 0.68  4
03 1 0 65 0.72 20
04 1 0 46 0.76 42
05 1 0 58 0.76  9
06 1 0 42 0.80 32
07 1 0 48 0.80 66
08 1 0 58 0.80 28
09 1 0 50 0.82 10
10 1 0 80 0.82  4
11 1 0 63 0.84 12
12 1 0 65 0.84 17
13 1 0 70 0.84 12
14 1 0 70 0.84 23
15 1 0 72 0.84 40
16 1 0 97 0.84 18
17 1 0 46 0.88 10
18 1 0 56 0.88 38
19 1 0 70 0.88  7
20 1 0 70 0.88 23
21 1 0 72 0.88 36
22 1 0 76 0.88  9
23 1 0 90 0.88 21
24 1 0 76 0.92 62
25 1 0 92 0.92 36
01 1 1 21 0.68 23
02 1 1 40 0.68 62
03 1 1 44 0.72 28
04 1 1 54 0.76 18
05 1 1 36 0.78 10
06 1 1 40 0.80 28
07 1 1 56 0.80 22
08 1 1 60 0.80 29
09 1 1 48 0.84 15
10 1 1 53 0.84 73
11 1 1 60 0.84 10
12 1 1 60 0.84  5
13 1 1 65 0.84 13
14 1 1 68 0.84 27
15 1 1 60 0.88 20
16 1 1 81 0.88 21
17 1 1 81 0.88 12
18 1 1 48 0.90 49
19 1 1 48 0.90 17
20 1 1 56 0.90 22
21 1 1 68 0.90 71
22 1 1 75 0.90 17
23 1 1 81 0.90 10
24 1 1 48 0.92 24
25 1 1 68 0.92 18
01 8 1 16 0.64 34
02 8 1 19 0.64  6
03 8 1 19 0.68  4
04 8 1 32 0.72 22
05 8 1 33 0.72 28
06 8 1 33 0.74 31
07 8 1 30 0.76 16
08 8 1 42 0.76 27
09 8 1 42 0.76  8
10 8 1 33 0.78 32
11 8 1 26 0.80 20
12 8 1 30 0.80 35
13 8 1 40 0.82 12
14 8 1 54 0.82 14
15 8 1 34 0.84 17
16 8 1 34 0.84 29
17 8 1 47 0.84 31
18 8 1 47 0.84  6
19 8 1 42 0.88 30
20 8 1 47 0.88 27
21 8 1 54 0.88 40
22 8 1 54 0.88 19
23 8 1 56 0.88  8
24 8 1 60 0.88  8
25 8 1 44 0.92 15
end


* create a event/censoring indicator for everyone

gen event = 1

* Quintiles

gen thorax_q = thorax

recode thorax_q 0/0.785=1 0.785/0.825= 2 0.825/0.845=3 0.845/0.885=4 0.885/0.99=5

* scenario where * staggered entry 0 to 42 days, funding runs out day 84 

gen day0 = mod(id,5)*10 + mod(id,3)
gen day1 = day0 + lngevity

gen t_fu  = lngevity*(day1 <= 84) + (84-day0)*(day1 > 84)
gen final = 1*(day1 <= 84)        +         0*(day1 > 84)

* summaries for 2 groups of interest

sort active

by active: summarize lngevity thorax sleep if partners==1  

* set up survival data

stset lngevity, failure(event)

*proportional hazards model

stcox active thorax if partners==1, nohr

stcox

*proportional hazards model, thorax-q values as strata 

stcox active if partners==1, strata(thorax_q) nohr

stcox


* ------------------------------------------

* set up survival data in 'staggered entry, funding runs out' scenario

stset t_fu, failure(final)

*proportional hazards model

stcox active thorax if partners==1, nohr

stcox


* -----------------------------------------

* Full followup, RA enters data days 6, 13, ..  makes a
* separate record for each week, or part therof

* unique id

gen uniqueid = 1000000 + partners*10000 + active*100 + id

keep uniqueid partners  active  lngevity  thorax  sleep  thorax_q

gen day_0  = mod(uniqueid,5) + 5*(mod(uniqueid,5)==0)
gen day_1  = day_0 + lngevity

expand 15   /* this makes 15 copies of the data, 1 per week */

sort uniqueid

gen e_day = 6 + 7 * ( mod(_n,15) + 15*(mod(_n,15)==0) -1 )

gen t0 = max(0,(e_day - day_0)-7)

gen t1 = min(e_day - day_0,lngevity)
    
gen complete = t1==lngevity

keep if (t0 < t1) & partners==1  /* this keeps the relevant weeks */

drop day_0 day_1 e_day

sort uniqueid t0
 
stset t1, failure(complete) enter(time t0) /* note syntax */

stcox active thorax if partners==1,  nohr

stcox

* --------

* the analysis by the boss

stcox active thorax if ( partners==1 & complete==1) ,  nohr

stcox