[R-sig-ME] Mixed model and repeated measures in R

[Thierry Onkelinx]

Dear Despina,

You have complete separation in your dataset. It shows in the output of the
GCA data. Extreme random intercept variances (SCAN_DATE:ID and ID), extreme
fixed effect parameters (intercept and Comb_PH_tod).

Your model is too complex for your data.

Best regards,

ir. Thierry Onkelinx
Statisticus / Statistician

Vlaamse Overheid / Government of Flanders
INSTITUUT VOOR NATUUR- EN BOSONDERZOEK / RESEARCH INSTITUTE FOR NATURE AND
FOREST
Team Biometrie & Kwaliteitszorg / Team Biometrics & Quality Assurance
thierry.onkelinx using inbo.be
Havenlaan 88 bus 73, 1000 Brussel
www.inbo.be

///////////////////////////////////////////////////////////////////////////////////////////
To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to say
what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner
The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey
///////////////////////////////////////////////////////////////////////////////////////////

<https://www.inbo.be>


Op zo 16 jun. 2019 om 21:13 schreef DESPINA MICHAILIDOU <
de.michailidou using gmail.com>:

> Hi All,
>
> I am trying to run regression analysis adjusted for repeated measures in R.
> The imaging pathology finding defined as Vert_effect, CA_effect,
> Vert_Intens etc is the outcome variable whereas the clinical symptom
> defined as Comb_PH_tod, Comb_PNP_tod etc, is the predictor variable. Other
> predictor variables that I am using are the daily prednisone use (Pred) and
> the use of immunosuppresive therapy (Immune_Categorical) or not. As the
> prednisone variable is being read as character in R i converted it to
> numeric because it is a number. For example some patients are getting 2 mg
> of prednisone but some others 40 mg.  I have two subset of diagnoses,  the
> one is TAK and the second one is GCA. As some patients have either right
> side posterior headache or left side posterior headache or both or none and
> either right side vertebral intensity (imaging study pathology) or left
> side vertebral intensity, or both or none vertebral intensity, for each
> patient I created two rows per subject. The first row represents the right
> sided symptoms and imaging pathology findings and the second row represents
> the left sided symptoms and imaging findings. My repeated measures are the
> side of the symptoms and imaging findings (had to create a separate
> variable for the right and left side symptoms and right and left imaging
> findings that I called it Side and put in there R, L, R, L etc), the ID of
> the patients and the Scan_date visit. Some patients had one scan visit but
> some other patients had multiple scan visits, and that is why I am
> considering scan visit date as a repeated measure. *So this is the code
> that i am using and for the subset of TAK i get this output*
>
> > TAK_data <- subset(Despina, Diagnosis=="TAK")
> > glmm_Vert_Intes <- glmer (Vert_Intes ~ Comb_PH_tod  + (1 |
> ID/SCAN_DATE/Side) + Pred + Immune_Catagorical , data=TAK_data,
> family=binomial(link = "logit"))
> Error in length(value <- as.numeric(value)) == 1L :
>   (maxstephalfit) PIRLS step-halvings failed to reduce deviance in
> pwrssUpdate
> > summary(glmm_Vert_Intes)
>
> *whereas for the subset of GCA patients there are no issues.*
> Generalized linear mixed model fit by maximum likelihood (Laplace
> Approximation) ['glmerMod']
>  Family: binomial  ( logit )
> Formula: Vert_Intes ~ Comb_PH_tod + (1 | ID/SCAN_DATE/Side) + Pred +
> Immune_Catagorical
>    Data: GCA_data
>
>      AIC      BIC   logLik deviance df.resid
>     87.8    113.0    -36.9     73.8      263
>
> Scaled residuals:
>      Min       1Q   Median       3Q      Max
> -0.98336 -0.00342 -0.00241 -0.00214  1.02148
>
> Random effects:
>  Groups              Name        Variance Std.Dev.
>  Side:(SCAN_DATE:ID) (Intercept)   0.00    0.00
>  SCAN_DATE:ID        (Intercept) 528.06   22.98
>  ID                  (Intercept)  18.84    4.34
> Number of obs: 270, groups:  Side:(SCAN_DATE:ID), 270; SCAN_DATE:ID, 135;
> ID, 54
>
> Fixed effects:
>                     Estimate Std. Error z value Pr(>|z|)
> (Intercept)        -10.63649    2.36017  -4.507 6.59e-06 ***
> Comb_PH_tod         -8.53678    4.56601  -1.870   0.0615 .
> Pred                -0.02353    0.09323  -0.252   0.8008
> Immune_Catagorical  -1.41376    2.69420  -0.525   0.5998
> ---
> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> Correlation of Fixed Effects:
>             (Intr) Cm_PH_ Pred
> Comb_PH_tod  0.299
> Pred        -0.289 -0.001
> Immn_Ctgrcl -0.520 -0.068  0.041
> convergence code: 0
> boundary (singular) fit: see ?isSingular
>
> *And this is the detailed code that I am using in R*
> install.packages("lme4")
> install.packages("readr")
>
> library(readr)
> library("lme4")
>
> setwd("~/Desktop/Despina")
> Despina <- read_csv("Despina.csv")
> as.factor(Despina$ID)
> as.factor(Despina$Diagnosis)
> as.factor(Despina$SCAN_DATE)
> as.factor(Despina$Side)
> as.factor(Despina$Immune_Catagorical)
> as.factor(Despina$LH_today)
> as.factor(Despina$PLH_today)
> as.factor(Despina$Dizz_today)
> as.factor(Despina$P_Diz_today)
> as.factor(Despina$CD_tod)
> as.factor(Despina$Head_today)
> as.factor(Despina$Vertig_today)
> as.factor(Despina$FTH_tod)
> as.factor(Despina$Comb_PH_tod)
> as.factor(Despina$Comb_NP_tod)
> as.factor(Despina$Comb_ANP_tod)
> as.factor(Despina$Comb_PNP_tod)
> as.factor(Despina$CNS_ever)
> as.factor(Despina$ULC_today)
> as.factor(Despina$Vert_effect)
> as.factor(Despina$CA_effect)
> as.factor(Despina$Sub_invol)
> as.factor(Despina$Ax_involv)
> as.factor(Despina$CA_intens)
> as.factor(Despina$Sub_intens)
> as.factor(Despina$Vert_Intes)
> as.factor(Despina$Ax_intens)
> as.factor(Despina$Comb_Vis_L_today)
> Despina$Pred<-as.numeric(as.character(Despina$Pred))
>
>
> TAK_data <- subset(Despina, Diagnosis=="TAK")
>
> glmm_Vert_Intes <- glmer (Vert_Intes ~ Comb_PH_tod  + (1 |
> ID/SCAN_DATE/Side) + Pred + Immune_Catagorical , data=TAK_data,
> family=binomial(link = "logit"))
> summary(glmm_Vert_Intes)
>
> GCA_data <- subset(Despina, Diagnosis=="GCA")
>
> glmm_Vert_Intes <- glmer (Vert_Intes ~ Comb_PH_tod  + (1 |
> ID/SCAN_DATE/Side) + Pred + Immune_Catagorical , data=GCA_data,
> family=binomial(link = "logit"))
> summary(glmm_Vert_Intes)
>
> *So my question is why i am getting this error in TAK patients and not in
> GCA patients?*
> Error in length(value <- as.numeric(value)) == 1L :
>   (maxstephalfit) PIRLS step-halvings failed to reduce deviance in
> pwrssUpdate
>
> Thank you all for your time and consideration in advance.
>
> Sincerely,
> Despina
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>

	[[alternative HTML version deleted]]