LME output changing overnight, can't restore

gm1 · September 2, 2023, 12:44pm

SOLVED

Hello! I'm working on my Msc dissertation and have been using linear mixed effect modeling in R for my project. As of 8/31 everything was working great, the Rstudio portion of my project was done, and I just had to go on to write up my analysis. As of 9/1 when I opened and ran the program to keep writing, the output entirely change. I don't know why, I've duplicated the file and gone through it doing any edits I can think of and nothing has changed it. I started from scratch and rewrote only the LME portion of the analysis. I even went through and re-cleaned my data file thinking maybe that got corrupted. My only hint is that a few lines of code that were unrelated to the output that's not working weren't there when I opened the file (a few lines plotting the density of the data). I've looked up documentation and seen that Rstudio has a file for the history of the file but I can't seem to figure out how to restore it? Or if I have, it's not doing anything. Additionally, I've reinstalled and reupdated the lme4 and lmerTest packages. My dissertation deadline is very close, please help!

relevant packages: lme4 and lmerTest

relevant code (I made it into a function, but the actual lme stuff was written per my supervisor's direction):

   myoptim <- "bobyqa"
   cd_opt <- FALSE

   control1 =glmerControl(optimizer=myoptim, optCtrl=list(maxfun=2e5), calc.derivs=cd_opt)

   InteractionModel <- glmer(DV ~ IV1 * IV2 + (1 | DF$P_ID) + (1 | face_stimulus), data = DF, control = control1,                             
                                                           family = inverse.gaussian(link = 'identity')

   NonInteractionModel <- glmer(DV ~ IV1 + IV2 + (1 | DF$P_ID) + (1 | face_stimulus), data = DF, 
                                                                  control = control1, family = inverse.gaussian(link = 'identity')

   header_func("Interaction")
   print(summary(InteractionModel))

   header_func("No Interaction")
   print(summary(NonInteractionModel))

   header_func("Chi-squared comparing models")
   anova(InteractionModel, NonInteractionModel)

All of this originally worked and had significant outputs. I've gone through several versions of this code and had the same output:

============================================================================"
[1] " "
[1] "Interaction"
[1] " "
[1] "----------------------------------------------------------------------------"
Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: inverse.gaussian  ( identity )
Formula: DV ~ IV1 * IV2 + (1 | DF$P_ID) + (1 | DF$face_stimulus)
   Data: DF
Control: control

     AIC      BIC   logLik deviance df.resid 
 14231.2  14273.2  -7106.6  14213.2      781 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.1650 -0.4817 -0.1722  0.3721  8.4469 

Random effects:
 Groups           Name        Variance  Std.Dev. 
 DF$P_ID          (Intercept) 2.931e+06 1.712e+03
 DF$face_stimulus (Intercept) 1.981e+04 1.407e+02
 Residual                     2.110e-04 1.453e-02
Number of obs: 790, groups:  DF$P_ID, 33; DF$face_stimulus, 24

Fixed effects:
            Estimate Std. Error t value Pr(>|z|)    
(Intercept)  3833.34     329.73  11.626  < 2e-16 ***
IV11         -519.12      79.03  -6.568 5.08e-11 ***
IV21          623.46     138.28   4.509 6.52e-06 ***
IV22         -189.10     112.01  -1.688   0.0914 .  
IV11:IV21     -20.91     131.62  -0.159   0.8738    
IV11:IV22    -182.57     104.45  -1.748   0.0805 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
          (Intr) IV11   IV21   IV22   IV11:IV21
IV11      -0.116                               
IV21       0.089 -0.210                        
IV22      -0.028 -0.027 -0.693                 
IV11:IV21 -0.051  0.484 -0.542  0.427          
IV11:IV22  0.001 -0.155  0.432 -0.687 -0.718   
[1] "============================================================================"
[1] " "
[1] "No Interaction"
[1] " "
[1] "----------------------------------------------------------------------------"
Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: inverse.gaussian  ( identity )
Formula: DV ~ IV1 + IV2 + (1 | DF$P_ID) + (1 | DF$face_stimulus)
   Data: DF
Control: control

     AIC      BIC   logLik deviance df.resid 
 14235.8  14268.5  -7110.9  14221.8      783 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.1477 -0.4852 -0.1755  0.3675  8.3921 

Random effects:
 Groups           Name        Variance  Std.Dev. 
 DF$P_ID          (Intercept) 2.901e+06 1.703e+03
 DF$face_stimulus (Intercept) 2.064e+04 1.437e+02
 Residual                     2.175e-04 1.475e-02
Number of obs: 790, groups:  DF$P_ID, 33; DF$face_stimulus, 24

Fixed effects:
            Estimate Std. Error t value Pr(>|z|)    
(Intercept)  3814.24     328.17  11.623  < 2e-16 ***
IV11         -492.57      67.71  -7.275 3.47e-13 ***
IV21          617.34     118.79   5.197 2.03e-07 ***
IV22         -343.15      82.17  -4.176 2.97e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
     (Intr) IV11   IV21  
IV11 -0.096              
IV21  0.077  0.073       
IV22 -0.057 -0.060 -0.711
[1] "============================================================================"
[1] " "
[1] "Chi-squared comparing models"
[1] " "
[1] "----------------------------------------------------------------------------"
Data: DF
Models:
NonInteractionModel: DV ~ IV1 + IV2 + (1 | DF$P_ID) + (1 | DF$face_stimulus)
InteractionModel: DV ~ IV1 * IV2 + (1 | DF$P_ID) + (1 | DF$face_stimulus)
                    npar   AIC   BIC  logLik deviance  Chisq Df Pr(>Chisq)  
NonInteractionModel    7 14236 14268 -7110.9    14222                       
InteractionModel       9 14231 14273 -7106.6    14213 8.5829  2    0.01369 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

But since opening it the other day, I keep getting stuck getting this output when I pass the same values to it:

[1] "============================================================================"
[1] " "
[1] "Interaction"
[1] " "
[1] "----------------------------------------------------------------------------"
Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: inverse.gaussian  ( identity )
Formula: DV ~ IV1 * IV2 + (1 | DF$P_ID) + (1 | DF$face_stimulus)
   Data: DF
Control: control1

     AIC      BIC   logLik deviance df.resid 
 14238.9  14280.9  -7110.4  14220.9      781 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.1235 -0.4676 -0.1599  0.3710  7.9442 

Random effects:
 Groups           Name        Variance  Std.Dev. 
 DF$P_ID          (Intercept) 1.554e+06 1.247e+03
 DF$face_stimulus (Intercept) 2.896e+04 1.702e+02
 Residual                     2.270e-04 1.507e-02
Number of obs: 790, groups:  DF$P_ID, 33; DF$face_stimulus, 24

Fixed effects:
            Estimate Std. Error t value Pr(>|z|)    
(Intercept)  3699.10     254.51  14.534  < 2e-16 ***
IV11         -565.10      85.52  -6.608 3.89e-11 ***
IV21          638.72     150.46   4.245 2.19e-05 ***
IV22         -186.13     123.91  -1.502   0.1331    
IV11:IV21     -11.47     141.67  -0.081   0.9355    
IV11:IV22    -201.66     113.83  -1.772   0.0765 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
          (Intr) IV11   IV21   IV22   IV11:IV21
IV11      -0.175                               
IV21       0.125 -0.202                        
IV22      -0.034 -0.037 -0.686                 
IV11:IV21 -0.071  0.469 -0.550  0.429          
IV11:IV22 -0.003 -0.135  0.437 -0.694 -0.715   
[1] "============================================================================"
[1] " "
[1] "No Interaction"
[1] " "
[1] "----------------------------------------------------------------------------"
Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: inverse.gaussian  ( identity )
Formula: DV ~ IV1 + IV2 + (1 | DF$P_ID) + (1 | DF$face_stimulus)
   Data: DF
Control: control

     AIC      BIC   logLik deviance df.resid 
 14235.8  14268.5  -7110.9  14221.8      783 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.1477 -0.4852 -0.1755  0.3675  8.3921 

Random effects:
 Groups           Name        Variance  Std.Dev. 
 DF$P_ID          (Intercept) 2.901e+06 1.703e+03
 DF$face_stimulus (Intercept) 2.064e+04 1.437e+02
 Residual                     2.175e-04 1.475e-02
Number of obs: 790, groups:  DF$P_ID, 33; DF$face_stimulus, 24

Fixed effects:
            Estimate Std. Error t value Pr(>|z|)    
(Intercept)  3814.24     328.17  11.623  < 2e-16 ***
IV11         -492.57      67.71  -7.275 3.47e-13 ***
IV21          617.33     118.79   5.197 2.03e-07 ***
IV22         -343.15      82.17  -4.176 2.97e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
     (Intr) IV11   IV21  
IV11 -0.096              
IV21  0.077  0.073       
IV22 -0.057 -0.060 -0.711
[1] "============================================================================"
[1] " "
[1] "Chi-squared comparing models"
[1] " "
[1] "----------------------------------------------------------------------------"
Data: DF
Models:
NonInteractionModel: DV ~ IV1 + IV2 + (1 | DF$P_ID) + (1 | DF$face_stimulus)
InteractionModel: DV ~ IV1 * IV2 + (1 | DF$P_ID) + (1 | DF$face_stimulus)
                    npar   AIC   BIC  logLik deviance  Chisq Df Pr(>Chisq)
NonInteractionModel    7 14236 14268 -7110.9    14222                     
InteractionModel       9 14239 14281 -7110.4    14221 0.8953  2     0.6391

What the heck?? Any help would be appreciated.

Edit: I should probably mention that this continued happening even when I used a version of the project from the 23rd with the current data.

technocrat · September 2, 2023, 9:47pm

The history file may be in .Rhistory in your working directory or in the same name under your home directory (I think it depends on whether the working directory has an .Rproj-user file (meaning an RStudio project had been created).

My experience with overnight mysteries like this

The reason it stopped working is that it was running in a stale session, so lots of stuff from earlier work was still in scope
Shoulda been backed up some way, such as git or a reprex

gm1 · September 3, 2023, 2:11pm

Hello! Thank you so much for your feedback. After a 48 hour process of Theseus' ship-ing it.. it turns out the issue was core files on my computer being corrupted and not R. Your reply definitely helped me narrow down what it could have been, but it turns out a forced wipe of the machine and fresh install of the OS was what was needed.

system · October 15, 2023, 2:11pm

This topic was automatically closed 42 days after the last reply. New replies are no longer allowed.

If you have a query related to it or one of the replies, start a new topic and refer back with a link.