[R-sig-ME] Fixed Effect and survival analysis (coxme package)
Sarahí Rueda Salazar
@rued@ @ending from ced@u@b@e@
Thu Jun 28 15:45:10 CEST 2018
Dear R users,
My question is related to the recently updated version of coxme package (
mixed effect in survival analysis). Currently, I work with multistate
models with changes in health status on elder population. I use coxPH by
stratifying the hazard by different transition type ( deterioration, health
improvements, death from healthy status and death from not healthy status).
So, I wondered if I can apply mixed effect on my multistate model to see
the differences by countries (shared frailty). I´ve read the coxme´s
document post on May 11, 2018 and others several materials in this subject
:Steele et al:2004, Austin : 2017,Putter: 2014, Willekens: 2014,
Brostön:2011, Mills, 2011 and finally Allison (which provide me a quite
understandable approach in the line of my humble knowledge in mathematical
demography) .
I thought it is not possible to use to multilevel in
multistate with coxme package so, I´ve done a model for each transition
type by separate . It means one model by transition type. If I stratify the
hazard (reference risk by interested covariates ) by my four transitions
type, I make hazard be free(no proportional by transition type) and it
means that the model estimates separately baseline hazard for the different
values of transition type ( following this material by Putter( 2018:7
<https://cran.r-project.org/web/packages/mstate/vignettes/Tutorial.pdf>)
e.g: It is a simple model using sex , reference category "male" with my
data:
> modelSex.0 <- coxph(Surv(Tstarta,Tstopa,status) ~
+ SexF.1+ SexF.2+SexF.3+SexF.4+
+ strata(trans),
+ data=d0)
SexF.1= related to covariate sex in transition 1(healthy to not healthy)
SexF.2 = transition healthy to death
SexF.3= transition Not healthy to Healthy
SexF.4= transition Not healthy to Death
But, I came across this post
<https://stat.ethz.ch/pipermail/r-help/2014-September/421690.html> where
Therneau replied a request ( I copy a chunk that I´m interested in):
2. The model above is the correct covariance structure for a set of
> families. There is a
> single intercept per subject, with a complex correlation matrix. The
> simpler "per family"
> frailty model would be
> model4 <- coxme(Surv(Survival, Event) ~ Sex + strata(cohort) + SNP1 + SNP2
> + SNP3 +
> (1|famid), death.dat)
> This model lets each family have a separate risk, with everyone in the
> same family sharing
> the exact same risk. It is less general than model3 above which lets a
> family have higher
> risk plus has variation between family members. A model with both
> per-subject and per family terms is identical to one with a covariance
> matrix of s1 K + s2 B, where K is the kinship matrix, B is a block
> diagonal matrix which
> has a solid block of "1" for each family, and s1 s2 are the fitted
> variance coefficients.
taking the example to my data, strata would be my transition type (trans)
instead "cohort" and my groups (instead famid) would be country, as follows
> modelSex.1 <- coxme(Surv(Tstarta,Tstopa,status) ~
+ SexF.1+ SexF.2+SexF.3+SexF.4+
+ (1|Country),
+ data=d0)
> anova(modelSex.0,modelSex.1)
Analysis of Deviance Table
Cox model: response is Surv(Tstarta, Tstopa, status)
Model 1: ~ SexF.1+ SexF.2+SexF.3+SexF.4+ strata(trans)
Model 2: ~ SexF.1+ SexF.2+SexF.3+SexF.4+ strata(trans) + (1 | Country)
loglik Chisq Df P(>|Chi|)
1 -1342082
2 -1339605 4953.8 1 < 2.2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> stem(exp(ranef(modelSex.1)[[1]]))
The decimal point is 1 digit(s) to the left of the |
6 | 589
8 | 1448477
10 | 115669
12 | 0846
> exp(ranef(modelSex.1)[[1]])
AT BE BG CY CZ EE
1.1981346 0.9712875 0.8092424 1.3614033 0.8382588 1.0117071
EL ES HU IE IT LT
0.7514988 1.2761863 1.0085404 0.8760514 1.1615717 0.9708070
LU LV MT PL PT RO
1.1893096 1.3381012 1.0522161 0.7843473 1.1642434 0.7911292
SK UK
0.9440166 0.8428222
> fixed.effects(modelSex.1)
SexF.1 SexF.2 SexF.3 SexF.4
0.16349517 -0.63184370 -0.08578787 -0.61023260
The thing is that I´m not sure on how to interpret the frailty values by
countries (random effect described by the variance within groups) because
I have four different effects (each related with my transition types). I
know that by using strata for my transition type is the same as I applied 4
different cox model (related for each type of transition). If I
apply separated model I would obtain frailty random effect by groups
(countries) regarding specific transition but, doing this model with the
strata ( modelSex.1 ) with my 4 transition type at once, I do not know
what those frailty values are telling me about the four different type of
hazard.
For other side, I have other question related to the recent article of
Mixed Effect (May, 2018) . Might be, it is a very very silly question but I
need to understand it .
In pag 9, second parraph , Therneau describes the simple cox model and it
is described the standard deviation (excess of risk for each group) . So,
there is an argument that states"... 15% of the families to be 1 std dev
or more above the mean..." Im working with that data and code but I
couldnt find where you got the value of 15%.
I hope somebody can help me,
Many thanks in advance,
Best wishes
Sarahí
*Centre d'Estudis DemogràficsCarrer de Ca n'Altayó, Edifici E2Universitat
Autònoma de Barcelona,08193
Bellaterra, Barcelona/SPAINPhone: 34/93.581.30.60
<34%2F93.581.30.60>Fax: 34/93.581.30.61
<34%2F93.581.30.61>e-mail: srueda using ced.uab.es
<ssancho using ced.uab.es>http://www.ced.uab.es <http://www.ced.uab.es/>*
--
*Sarahí Rueda Salazar*
*Investigadora en Formación (FPI/CED)*
*Centre d'Estudis DemogràficsCarrer de Ca n'Altayó, Edifici E2Universitat
Autònoma de Barcelona,08193
Bellaterra, Barcelona/SPAINPhone: 34/93.581.30.60
<34%2F93.581.30.60>Fax: 34/93.581.30.61
<34%2F93.581.30.61>e-mail: srueda using ced.uab.es
<ssancho using ced.uab.es>http://www.ced.uab.es <http://www.ced.uab.es/>*
More information about the R-sig-mixed-models
mailing list