Fine and Gray with Censoring

[Caecilia Lipowski]

When I produce cumulative incidence curves (CIC) for the Fine and Gray model using stcrreg and stcurve, I noticed that the estimated plateau level of the CIC (for large durations) highly depends on the degree of censoring – see the graph below. E.g. if 60% of the spells are censored, the CIC attains a maximum of 0.18, if instead 40% are censored, the maximum CIC is 0.29. Censoring is random. In theory, however, the estimated plateau level of the CIC should NOT depend on the degree of censoring. Indeed, the estimated plateau level of the CIC is invariant, when I do the same exercise in R using "cmprsk". So, either, I do a basic coding error or there is an issue with the Stata implementation.

Here is the Stata code generating the data plus the graph:

set seed 1000

foreach level in 40 60 80 100 { // level of censoring

clear
set obs 1000
gen i=_n
gen ra=round(runiform())
gen rb=1-ra
gen udur=-ln(uniform())

gen female=rnormal(0,1)
gen deutsch=rnormal(0,1)
gen VT_HE=rnormal(0,1)

*add censoring to data:
qui gen ran=runiform()
qui gen ran2=runiform()
qui replace ra=0 if ran>0.`level'
qui replace rb=0 if ran>0.`level'
qui replace udur=udur*ran2 if ran>0.`level'

stset udur, failure(ra==1)
stcrreg VT_HE deutsch female, compete(ra==0) iter(100)

stcurve, cif legend(off) xla(, grid) outfile("FG_`level'.dta", replace)
}

foreach level in 40 60 80 100 { // level of censoring
use "FG_`level'.dta", clear
rename ci1 ci`level'
save "FG_`level'.dta", replace
}

use "FG_40.dta", clear
foreach level in 60 80 100 {
joinby _t using FG_`level'.dta, unmatched(both)
drop _merge
}

sort _t
twoway (line ci40 _t)(line ci60 _t)(line ci80 _t)(line ci100 _t), legend(order(1 "60" 2 "40" 3 "20" 4 "0"))
graph export statalist.png, replace


foreach level in 40 60 80 100 {
erase "FG_`level'.dta"
}

Has anybody come across this problem before?
Thank you in advance!

[Steve Samuels]

Hello, Caecelia!

Thanks for the nice graph. It illustrated your problem better than any words.

What you are seeing is expected behavior, but before I explain, let me clear up some terminology.
Your use of the word “censoring” is not correct. In a two-cause competing risk analysis, there are three types of outomes.

1. The person died of the cause of interest.
2. The person died of from some other cause (the competing risk).
3. The person was alive at the end of observation. This and only this outcome is a censored observation.

You have conflated types 2 and 3. Unfortunately, your code reflects this confusion. stset designates a death from the cause of interest with ra = 1. To Stata, that means that ra = 0 is a censored observation. However stcrreg designates death from a competing cause with the same indicator: ra = 0. Perhaps you intended one of ra and rb to be the censoring indicator. This would have led to the same problem: ra and rb are identical.

Your question:

The phenomenon you observe is expected. In the Fine-Gray equations, observations with the competing event remain in the risk set. Therefore, as the number of competing events increases, the number at risk for death from the cause of interest increases. This quantity, very roughly, is the denominator for the CIF. As you hold constant the numerator (deaths from cause of interest), the CIF decreases.. A good review is Austin et al, 2016.

P.S. I found your code somewhat hard to follow. You created terms rb and ran2, but never used them. In the future, adding comments and labeling variables would be helpful.


References:

Jason P. Fine and Robert J. Gray, A Proportional Hazards Model for the Subdistribution of a Competing Risk
Journal of the American Statistical Association, Vol. 94, No. 446 (Jun. 1999), pp. 496- 509

Austin, P. C., Lee, D. S.and Fine, J. P. (2016). Introduction to the Analysis of Survival Data in the Presence of Competing Risks. Circulation, 133(6), 601-609

[Steve Samuels]

I must apologize,Caecelia. You certainly know what censoring means. The problem is the error in your stcrreg statment.:
compete(ra==0) where ra = 0 also is your censoring indicator. So your foreach code not only created differing censoring proportions, but also matching proportions of observations with competing events

[Caecilia Lipowski]

Hi Steve,

Thank you for solving this problem, your first reply was already very helpful! Changing "compete(ra==0)" to "compete(rb==1)" did the job by untangling censored from competing events.