Hi everyone,

Longtime Statelist forum lurker and clinician here who is a PhD student with an interest in clinical trials but lacks competency with statistical software.
I’m hoping I can get some help with the survsim command developed by Michael J. Crowther (https://www.mjcrowther.co.uk).

I am using survsim in Stata 14.2 on Mac OS High Sierra 10.13.6

I am trying to do a simulation study to see the effect of run-in period duration and adherence cutoffs to exclude non-adherent participants on trial power.
I’m varying a couple of parameters but am working on the baseline case right now and won't go into detail beyond the baseline case.

I have been successful in generating the data and running the simulation but after reinstalling survsim I am getting an error which states “invalid syntax r(198);” but it actually still generates a survival time in a now column.

This had previously worked without any issues but in an attempt to add time-varying covariates to the model I downloaded the following packages: merlin, moremata and installed the updated simsurv 4.0 program (https://www.mjcrowther.co.uk/code/survsim) and uninstalled the older version that I downloaded from a different server.

Here’s some information about my data:
The baseline case has:
a) n=2750
b) proportion of total non-adherers (i.e. 0 study drugs taken) = 0.05
c) an 8 week run-in period with baseline adherence = to 90% that decreases 1% per week
(l generated longitudinal data with a random slope of -1/week and intercept 90%)
d) adherence cutoff of 80% with dichotimization into adherence/non-adherent
e) proportion of hyperkalemia during run-in period leading to exclusion = 0.20
a run-in failure is defined as any participant with adherence<80% or hyperkalemia
f) participants that are randomized are allocated 1:1 to active or placebo
g) post-randomization adherence time non-varying and is equal to the run-in adherence at at 8 weeks which is used to generate survival times.
It is assumed that patients in the placebo have 0% adherence (they are not taking active study drug)…I did this because I haven’t been able to figure out interactions with survsim yet
The HR is 0.75 for the effect size but since it is modelled off of adherence instead which ranged from 0-100% so it is HR=0.0075 instead (linearity assumption)
I am also planning on changing adherence to be a time-dependent variable at some point...

Here is my data:
(I’ve only reported the 8 week values)… missing means a participant wasn’t randomized and was excluded because of non-adherence or hyperkalemia

input float(id runin_adherence8 runin_nonadherence hyperkalemia runin_failure8 randomization8 randomization8_treatment postrandomization_adherence8) double stime8
1 89.38921 0 1 1 0 . . .
2 74.60295 0 0 1 0 . . .
3 72.91611 0 0 1 0 . . .
4 97.27155 0 0 0 1 1 97.27155 .3169711239803543
5 90.70451 0 0 0 1 1 90.70451 21.20386339416367

This line of code that is specifically getting the error message is:
survsim stime8, distribution(exponential) lambdas(0.1) covariates (postrandomization_adherence8 -0.0075)

When I try other examples provided in ado files they all seem to work.

Does anybody have any idea why it is no longer working? Did something change in the syntax between versions that I am missing?
I don’t see any obvious syntax error in my code and unfortunately my simulation program won’t run with this error despite it generating survival times.

Thanks,

The full code is below:

*program
program baselinecase, rclass
*clear data
clear
*set the number of observations
set obs 2750
*generate id
gen id = _n
*runin non-adherers
gen runin_nonadherence=rbinomial(1, 0.05)
*longitudinal data for run-in adherence with random slope with intercept = 90%, B for visit = -1
generate u_0i = rnormal(0,3)
generate u_1i = rnormal(0,1)
expand 8
bysort id: generate runin_visit = _n
generate e_ij = rnormal(0,2)
generate runin_adherence = 90 -1*runin_visit + u_0i + runin_visit*u_1i + e_ij
drop u_0i u_1i e_ij
*reshaping data
reshape wide runin_adherence, i(id) j(runin_visit)
*changing runin_adherence in nonadherers
local runin_adherence runin_adherence1 runin_adherence2 runin_adherence3 runin_adherence4 runin_adherence5 runin_adherence6 runin_adherence7 runin_adherence8
foreach var of local runin_adherence{
replace `var'=0 if runin_nonadherence==1
}
*dichotomization
local runin_adherence runin_adherence1 runin_adherence2 runin_adherence3 runin_adherence4 runin_adherence5 runin_adherence6 runin_adherence7 runin_adherence8
foreach var of local runin_adherence{
gen `var'_binary=1 if `var'>=80
replace `var'_binary=0 if `var'<80
}
gen hyperkalemia=rbinomial(1, 0.20) if runin_nonadherence==0
replace hyperkalemia=0 if runin_nonadherence==1
*run-in failure
local runin_failure runin_failure1 runin_failure2 runin_failure3 runin_failure4 runin_failure5 runin_failure6 runin_failure7 runin_failure8
foreach var of local runin_failure{
gen `var'=.
}
replace runin_failure1=1 if runin_adherence1_binary==0 | hyperkalemia==1
replace runin_failure1=0 if runin_adherence1_binary==1 & hyperkalemia==0
replace runin_failure2=1 if runin_adherence2_binary==0 | hyperkalemia==1
replace runin_failure2=0 if runin_adherence2_binary==1 & hyperkalemia==0
replace runin_failure3=1 if runin_adherence3_binary==0 | hyperkalemia==1
replace runin_failure3=0 if runin_adherence3_binary==1 & hyperkalemia==0
replace runin_failure4=1 if runin_adherence4_binary==0 | hyperkalemia==1
replace runin_failure4=0 if runin_adherence4_binary==1 & hyperkalemia==0
replace runin_failure5=1 if runin_adherence5_binary==0 | hyperkalemia==1
replace runin_failure5=0 if runin_adherence5_binary==1 & hyperkalemia==0
replace runin_failure6=1 if runin_adherence6_binary==0 | hyperkalemia==1
replace runin_failure6=0 if runin_adherence6_binary==1 & hyperkalemia==0
replace runin_failure7=1 if runin_adherence7_binary==0 | hyperkalemia==1
replace runin_failure7=0 if runin_adherence7_binary==1 & hyperkalemia==0
replace runin_failure8=1 if runin_adherence8_binary==0 | hyperkalemia==1
replace runin_failure8=0 if runin_adherence8_binary==1 & hyperkalemia==0
*randomization
local randomization randomization1 randomization2 randomization3 randomization4 randomization5 randomization6 randomization7 randomization8
foreach var of local randomization{
gen `var'=.
}
replace randomization1=1 if runin_failure1==0
replace randomization1=0 if runin_failure1==1
replace randomization2=1 if runin_failure2==0
replace randomization2=0 if runin_failure2==1
replace randomization3=1 if runin_failure3==0
replace randomization3=0 if runin_failure3==1
replace randomization4=1 if runin_failure4==0
replace randomization4=0 if runin_failure4==1
replace randomization5=1 if runin_failure5==0
replace randomization5=0 if runin_failure5==1
replace randomization6=1 if runin_failure6==0
replace randomization6=0 if runin_failure6==1
replace randomization7=1 if runin_failure7==0
replace randomization7=0 if runin_failure7==1
replace randomization8=1 if runin_failure8==0
replace randomization8=0 if runin_failure8==1
*randomization to treatment groups
local randomization randomization1 randomization2 randomization3 randomization4 randomization5 randomization6 randomization7 randomization8
foreach var of local randomization{
gen `var'_treatment = runiform()>0.5 if `var'==1
}
*postrandomization adherence
gen postrandomization_adherence8 = runin_adherence8 if randomization8_treatment==1
replace postrandomization_adherence8=0 if randomization8_treatment==0
replace postrandomization_adherence8=100 if postrandomization_adherence8>100 & postrandomization_adherence8!=.
*simulating survival time
survsim stime8, distribution(exponential) lambdas(0.1) covariates (postrandomization_adherence8 -0.0075)
*events
generate event8=1 if stime8 <= 5
replace event8=0 if stime8>5
*right censoring
replace stime8 = 5 if event8 == 0
*Cox PH model,
stset stime8, failure(event8 = 1)
streg randomization8_treatment, distribution(exponential) nohr
return scalar hr = _b[randomization8_treatment]
return scalar sehr = _se[randomization8_treatment]
return scalar t = (_b[randomization8_treatment]/_se[randomization8_treatment])
return scalar p = (2 * ttail(e(df_r), abs(_b[randomization8_treatment]/_se[randomization8_treatment])))
end

simulate hr = r(hr) sehr = r(sehr) t = r(t) p = r(p), reps(1) seed(0) nodots nolegend: baselinecase
simulate hr = r(hr) sehr = r(sehr) t = r(t) p = r(p), reps(10) seed(0) nodots nolegend: baselinecase
simulate hr = r(hr) sehr = r(sehr) t = r(t) p = r(p), reps(100) seed(0) nodots nolegend: baselinecase
simulate hr = r(hr) sehr = r(sehr) t = r(t) p = r(p), reps(1000) seed(0) nodots nolegend: baselinecase
simulate hr = r(hr) sehr = r(sehr) t = r(t) p = r(p), reps(10000) seed(0) nodots nolegend: baselinecase