[Sean Hardiman]

Hi Clyde, sorry about the confusion. I should have included mention of study groups up front, my apologies! I'm still very much a novice with programming, are there any particular resources I might look at to figure out how to implement your suggestion of writing loops?

[Clyde Schechter]

You should familiarize yourself with the -foreach- and -forvalues- commands to learn about loops. The -levelsof- command is also important for looping over values of a variable (which is what you need to do here.) First open the help file for each of these commands and then click on the blue link near the top to get to the corresponding chapter in the PDF documentation. They are all well explained there.

That said, this problem is not a good one to use as your first application of loops. So here's a quick hack that will work. There are better ways to do it, but this is the simplest, requires the fewest changes to the code you hav, and will least burden you with learning unfamiliar aspects of Stata while you are working on a production project.

capture program drop get_survival
program define get_survival, rclass
    tempfile life_table
    sts list, risktable(30 365 730 1095) saving(`life_table')
    preserve
    use `life_table', clear
    return scalar S30 = survivor[1]
    return scalar S365 = survivor[2]
    return scalar S730 = survivor[3]
    return scalar S1095 = survivor[4]
    restore
    exit
end

// DON'T FORGET TO -STSET- YOUR DATA.
// AND BE SURE TO SET THE RANDOM NUMBER SEED TO GET REPRODUCIBLE RESULTS
levelsof study_group, local(groups)
foreach g of local groups {
    preserve
    keep if study_group == `g'
    bootstrap S30 = r(S30) S365 = r(365) S730 = r(S730) S1095 = r(S1095), reps(1000): ///
        get_survival
    estat bootstrap
    restore
}

[Jennifer Yo]

Dear Clyde,
This was discussion was incredibly helpful to my project where I had the exact same issue as Sean.
In addition to the survival, I need to obtain the absolute risk difference for different timepoints (eg. 3, 5, 8 years) with 95% confidence interval.
I'm very new to stata programming and would be interested in your suggestion of how to implement a code to achieve this?
Thank you

[Clyde Schechter]

Assuming your data set is structured in the same way as O.P.'s, it's a minor modification of the code:

capture program drop get_survival
program define get_survival, rclass
    tempfile life_table
    sts list, risktable(3 5 8) saving(`life_table')
    preserve
    use `life_table', clear
    return scalar S3 = survivor[1]
    return scalar S5 = survivor[2]
    return scalar S8 = survivor[3]
    return scalar D3_5 = survivor[1] - survivor[2]
    return scalar D3_8 = survivor[1] - survivor[3]
    return scalar D5_8 = survivor[2] - survivor[3]
    restore
    exit
end

// DON'T FORGET TO -STSET- YOUR DATA.
// AND BE SURE TO SET THE RANDOM NUMBER SEED TO GET REPRODUCIBLE RESULTS
levelsof study_group, local(groups)
foreach g of local groups {
    preserve
    keep if study_group == `g'
    bootstrap S3 = r(S3) S5 = r(S5) S8 = r(S8) D3_5 = r(D3_5) D3_8 = r(D3_8) ///
        D5_8 = r(D5_8), reps(1000): // GET SURVIVAL
    estat bootstrap
    restore
}

[Jennifer Yo]

Thank you Clyde!
Just to clarify what I meant to say was I need the absolute risk difference (between groups eg. exposed vs unexposed) at different time points (eg. 3,5, 8 years).
Apologies for the confusion.

[Clyde Schechter]

capture program drop get_survival
program define get_survival, rclass
    tempfile life_table
    sts list, risktable(3 5 8) by(group) saving(`life_table')
    preserve
    use `life_table', clear
    return scalar D3 = survivor[4] - survivor[1]
    return scalar D5 = survivor[5] - survivor[2]
    return scalar D8 = survivor[6] - survivor[3]
    restore
    exit
end

// DON'T FORGET TO -STSET- YOUR DATA.
// AND BE SURE TO SET THE RANDOM NUMBER SEED TO GET REPRODUCIBLE RESULTS
bootstrap D3 = r(D3) D5 = r(D5) D8 = r(D8), reps(1000): // GET SURVIVAL
estat bootstrap

[Jennifer Yo]

Thank you very much Clyde.
When I used your above code to obtain risk difference, I don't believe the by(group) in the sts line worked.
I am trying to obtain the absolute risk by groups (with bootstrapped CIs) and absolute risk difference between groups (with bootstrapped CIs).
Below is the code I used. I would greatly appreciate any insights as to where I have gone wrong in my code.

Thank you again,
Jennifer

*write bootstrapping program*
capture program drop get_absoluterisk
program define get_absoluterisk, rclass
tempfile life_table
sts list, failure risktable(36 72 108 144) saving (`life_table')
preserve
use `life_table', clear
return scalar f36= failure[1]
return scalar f72 = failure[2]
return scalar f108= failure[3]
return scalar f144= failure[4]
restore
exit
end


*estimate boostrap CI for absolute risk in pe expsoed and unexposed*
levelsof pe, local(groups)
foreach g of local groups {
preserve
keep if pe==`g'
bootstrap f36=r(f36) f72=r(f72) f108=r(f108) f144=r(f144), reps(2000) seed(12345): get_absoluterisk
restore
}


*obtain absolute risk difference*
capture program drop get_absoluteriskdiff
program define get_absoluteriskdiff, rclass
tempfile life_table
sts list, failure by(pe) risktable(36 72 108 144) saving (`life_table')
preserve
use `life_table', clear
return scalar D36= failure[5] - failure[1]
return scalar D72 = failure[6] - failure[2]
return scalar D108= failure[7] - failure[3]
return scalar D144= failure[8] - failure[4]
restore
exit
end


bootstrap D36=r(D36) D72=r(D72) D108=r(D108) D144=r(D144), reps(1000) seed(12345): get_absoluteriskdiff