Issue with stcurve in loop after estimates store/save

[Bree Shi]

Dear all,

I am analying a dataset on cancer using Fine and Gray model (stcrreg) and is trying to draw CIF on various cancer types after stcrreg.
Since it is taking a long time for stcrreg to run, I have been tryig to use estimates store/save on cancer type 1 to 13, then restore the estimates and then draw the curves over loop so I can modify how the graphs look without having to re-run the models. However, I have bumped into the issue of stcurve returning "no obseravtions" error after estimates restore/use. The code I am using are as follows.

foreach cat in  1 2 3 4 5 6 7 8 9 10 11 12 13 { 
    use "/Users/breeshi/Desktop/Cancer/Saved dataset/Analysis SET B_step3 Ready for analysis.dta", clear
    stset bmi_to, failure(cancercat == `cat') id(id) origin(bmi_from) scale(365.25)
    replace failuretype=5 if cancercat~=. & cancercat~=`cat'
    stcrreg i.obese_vary , compete(failuretype==2 3 5)
    local title  `:label (cancercat) `cat''
    estimates store cancercat`cat'
    estimates save `"cancercat`cat'"',replace
    }     

foreach cat in  1 2 3 4 5 6 7 8 9 10 11 12 13 { 
        estimates restore cancercat`cat'
        stcurve, at(obese_vary=(0 1)) cif ytitle("Cumulative Cancer Incidence",size(small)) ///
        graphregion(color(white)) title("`title'",size(med)) legend(off) ///
        xtitle("Years post-transplant",size(small)) lpattern( dash solid) ///
        lcolor(ebblue orange_red) xlabel(, labsize(small)) ylabel(, labsize(small)) saving("SetB_cancercat_`cat'", replace)
      }

The code gives me the first graph (i.e. cancercat==1) graph alright, but from the 2nd graph, it gives me "no observations
r(2000);" It is highly inefficient for me and I have been pulling my hair out trying to resolve this.

Appreciate any help in advance.

Best Regards,
Bree

[Clyde Schechter]

You don't show example data, but I tried doing a simulation of your problem using the StataCorp website hypoxia data set and a simplified version of your code. In doing that, I did not replicate the failure mode you describe. However, I will point out that in order for this to work properly, you must re-stset- your data for each graph so that it designates the appropriate failure types. The following works:

clear*
webuse hypoxia
levelsof failtype, local(failure_types)
foreach i of local failure_types {
    local failtype: list failure_types - i
    stset dftime, failure(failtype == `i')
    stcrreg ifp tumsize pelnode, compete(failtype == `failtype')
    estimates store fail`i'
    stcurve, cif
}

foreach i of local failure_types {
    local failtype: list failure_types - i
    stset dftime, failure(failtype == `i') // CODE WILL NOT RUN PROPERLY WITHOUT THIS
    estimates restore fail`i'
    stcurve, cif name(fail`i', replace)
}

Perhaps if you revise your code modeled on this approach you will have success. If not, when posting back, please show example data, using the -dataex- command.

If you are running version 18, 17, 16 or a fully updated version 15.1 or 14.2, -dataex- is already part of your official Stata installation. If not, run -ssc install dataex- to get it. Either way, run -help dataex- to read the simple instructions for using it. -dataex- will save you time; it is easier and quicker than typing out tables. It includes complete information about aspects of the data that are often critical to answering your question but cannot be seen from tabular displays or screenshots. It also makes it possible for those who want to help you to create a faithful representation of your example to try out their code, which in turn makes it more likely that their answer will actually work in your data.

[Bree Shi]

Many thanks for your reply, Clyde.

I tried your suggested code, but I am still encountering the "no observations r(2000)" error. Since the stcrreg model takes over 24 hours to run, it’s highly frustrating to have to rerun the entire model just to modify the graph appearance.

To help troubleshoot this issue, I have attached a simplified version of my dataset below. Thank you again for your assistance.

clear
input str9 id float(failuretype obese_vary bmi_from bmi_to) byte(_st _d) int _origin double(_t _t0)
"P1000463" 2 0 15194 17057 1 0 15194   5.100616016427105 0
"P1000518" 1 0 15448 19698 1 1 15448   11.63586584531143 0
"P1000530" 3 0 17049 21034 1 0 17049  10.910335386721425 0
"P1000631" 2 0 15480 21746 1 0 15480  17.155373032169745 0
"P1000788" 2 1 16102 20713 1 0 16102  12.624229979466119 0
"P1000811" 1 0 15263 18117 1 1 15263   7.813826146475017 0
"P1000890" 3 1 15291 16852 1 0 15291    4.27378507871321 0
"P1001092" 3 0 15270 20720 1 0 15270  14.921286789869953 0
"P1001148" 2 0 15127 15148 1 0 15127 .057494866529774126 0
"P1001159" 0 0 15850 23375 1 0 15850  20.602327173169062 0
end
format %td bmi_from
format %td bmi_to
label values failuretype failuretype
label def failuretype 0 "Censored", modify
label def failuretype 1 "Cancer incidence", modify
label def failuretype 2 "Graft failure", modify
label def failuretype 3 "No-cancer death", modify
label values obese_vary obese_vary
label def obese_vary 0 "Not obese", modify
label def obese_vary 1 "Obese", modify

The code I used is

    use "/Users/breeshi/Desktop/Cancer/Saved dataset/Analysis set B_step3 Ready for analysis.dta", clear
    keep in 1/1000
    stset bmi_to, failure(failuretype==1) id(id) origin(bmi_from) scale(365.25)
    stcrreg     i.obese_vary ,  compete(failuretype==2 3)        
    estimate save "/Users/breeshi/Desktop/Cancer/Saved dataset/test",replace

    
    use "/Users/breeshi/Desktop/Cancer/Saved dataset/Analysis set B_step3 Ready for analysis.dta", clear
    keep in 1/1000
    stset bmi_to, failure(failuretype==1) id(id) origin(bmi_from) scale(365.25)
    estimate use "/Users/breeshi/Desktop/Cancer/Saved dataset/test"
    stcurve, at(obese_vary=(0 1)) cif

.         use "/Users/breeshi/Desktop/Cancer/Saved dataset/Analysis set B_step3 Ready for an
> alysis.dta", clear

.         keep in 1/1000
(16,036 observations deleted)

.         stset bmi_to, failure(failuretype==1) id(id) origin(bmi_from) scale(365.25)

Survival-time data settings

           ID variable: id
         Failure event: failuretype==1
Observed time interval: (bmi_to[_n-1], bmi_to]
     Exit on or before: failure
     Time for analysis: (time-origin)/365.25
                Origin: time bmi_from

--------------------------------------------------------------------------
      1,000  total observations
          0  exclusions
--------------------------------------------------------------------------
      1,000  observations remaining, representing
        950  subjects
         44  failures in single-failure-per-subject data
  3,768.211  total analysis time at risk and under observation
                                                At risk from t =         0
                                     Earliest observed entry t =         0
                                          Last observed exit t =  23.72348


.     estimate restore test
(results test are active now)

.     stcurve, at(obese_vary=(0 1)) cif
no observations
r(2000);

I have also tried estimate store/restore and having the same error message.

[Clyde Schechter]

Ah, I think I have found it.

Whenever you run -estimates use-, Stata has a problem. You are bringing in a set of estimates, but Stata has no way to know what data were used to create those estimates. It has no way to know which observations are in the estimation sample and which are not. So, by default, it automatically sets the estimation sample to be empty. You have to use the -estimates esample:- command after -estimates use- runs to tell Stata what the estimation sample is.

    use "/Users/breeshi/Desktop/Cancer/Saved dataset/Analysis set B_step3 Ready for analysis.dta", clear
    keep in 1/1000
    stset bmi_to, failure(failuretype==1) id(id) origin(bmi_from) scale(365.25)
    estimate use "/Users/breeshi/Desktop/Cancer/Saved dataset/test"
    estimates esample:
    stcurve, at(obese_vary=(0 1)) cif

will do the trick. As specified, with no varlist, or if, or in conditions this tells Stata to use the entire data set in memory as the estimation sample. But you can specify any subset of the variables and any subset of the observations by appropriately writing a variable list, and -if- or -in- clause (or both) after the colon.

[Bree Shi]

Hurrah! This works. Thank you so much, Clyde.