1
$\begingroup$

I am trying to run the following model:

            I(week^2):mutation_status +
            (week + I(week^2) | subject_id) , 
          data = sim_dat)

This is the output I get from this model. There is one random slope that is rather high (-0.95) and I was curious to know how I can change my data simulation code to lower that value. I am also having issues with convergence and would love to get thoughts on how to avoid issues regarding convergence. Thank you!

Linear mixed model fit by REML. t-tests use Satterthwaite's method [lmerModLmerTest]

Formula: fetal_weight ~ week + mutation_status + week:mutation_status + I(week^2) +
 I(week^2):mutation_status + (week + I(week^2) | subject_id)
Data: sim_dat

REML criterion at convergence: 114

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.9343 -0.3980 -0.0075  0.4168  1.9928 

Random effects:
 Groups     Name        Variance  Std.Dev. Corr       
 subject_id (Intercept)   1.51524  1.2309             
            week         61.01832  7.8114   0.16      
            I(week^2)   496.20236 22.2756  -0.19 -0.95
 Residual                 0.05892  0.2427             
Number of obs: 100, groups:  subject_id, 20

Fixed effects:
                           Estimate Std. Error       df t value       Pr(>|t|)    
(Intercept)                  3.1331     0.4226  17.9961   7.413 0.000000713699 ***
week                        -1.2465     3.5210  18.0053  -0.354          0.727    
mutation_statusY             0.5061     0.5977  17.9961   0.847          0.408    
I(week^2)                   34.5706    10.8146  18.0084   3.197          0.005 ** 
week:mutation_statusY       -1.0697     4.9795  18.0053  -0.215          0.832    
mutation_statusY:I(week^2) 202.5521    15.2942  18.0084  13.244 0.000000000101 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
            (Intr) week   mttn_Y I(w^2) wk:m_Y
week        -0.160                            
muttn_sttsY -0.707  0.113                     
I(week^2)    0.145 -0.966 -0.103              
wk:mttn_stY  0.113 -0.707 -0.160  0.683       
mtt_Y:I(^2) -0.103  0.683  0.145 -0.707 -0.966

Code for data simulation:

set.seed(123)  # set the seed
J <- 20  # number of individuals (clusters)
cs <- 5  # number of time points (cluster size)
gam <- c(0, 0.75, 1.34)  # fixed effects
G <- matrix(c(1.75, 0, 0,
              0, 0.0125, 0,
              0,0,0.0625), 
            nrow = 3)  # random effect variances (G-matrix)
sigma2 <- 0.08

mutation_status <- rep(0:1, each = 5)
mutation_status <- rep(mutation_status, times = 10)

X <- cbind(1, seq_len(cs), 
           (seq_len(cs))^2)  # for each individual
X <- X[rep(seq_len(cs), J), ] 
X <- X[rep(seq_len(cs), J), ] 
X[,2] <- X[,2]*0.05
X[,3] <- X[,3]*0.05

# repeat each row cs times
pid <- seq_len(J)  # individual id
pid <- rep(pid, each = cs)


# Generate person-level (lv-2) random effects
uj <- lmf::rmnorm(J, mean = rep(0, 2), varcov = G)

# Generate repeated-measure-level (lv-1) error term
eij <- rnorm(J * cs, sd = sqrt(sigma2))

# Compute beta_j's
betaj <- matrix(gam, nrow = J, ncol = 3, byrow = TRUE) + uj

# Compute outcome:
y <- rowSums(X * betaj[pid, ]) + eij + mutation_status*X[ , 3]*10 + 3

# Form a data frame
sim_dat <- tibble(y, time = X[ , 2], pid, mutation_status)
sim_dat <- sim_dat[,c(3,2,4,1)]
colnames(sim_dat) <- c("subject_id", "week", "mutation_status", "fetal_weight")
sim_dat$week_2 <-sim_dat$week*5 + 15

sim_dat$mutation_status <- as.factor(ifelse(sim_dat$mutation_status == 0, "N", "Y"))
sim_dat$subject_id <- as.factor(sim_dat$subject_id )
#Model Building 
$\endgroup$
2
  • $\begingroup$ Convergence is a problem with your model, not a problem with the programming. Agree that this question is better suited to Cross Validated. $\endgroup$
    – Michelle
    Jul 10 '20 at 10:34
  • $\begingroup$ Apologies, I actually misspoke. I do not want the correlation of the random slopes between week and week^2 to be so highly correlated. But I will try to ask the question stackexchange instead! Thanks! $\endgroup$ Jul 10 '20 at 15:16

Browse other questions tagged or ask your own question.