Mixed model and repeated measures in R

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 at 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 at 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]]