[R-sig-ME] Fitted values in lmer

Henric Nilsson (Public) nilsson.henric at gmail.com
Fri Mar 14 22:28:37 CET 2008


Douglas Bates wrote:

> What I was referring to is the na.action attribute of the model frame
> 
>> library(lme4)
> Loading required package: Matrix
> Loading required package: lattice
> 
> Attaching package: 'Matrix'
> 
> 
> 	The following object(s) are masked from package:stats :
> 
> 	 xtabs
> 
>> slps <- sleepstudy
>> slps$Reaction[c(14, 34, 92)] <- NA
>> fm1 <- lmer(Reaction ~ Days + (Days|Subject), slps)
>> str(fm1 at frame)
> 'data.frame':	177 obs. of  3 variables:
>  $ Reaction: num  250 259 251 321 357 ...
>  $ Days    : num  0 1 2 3 4 5 6 7 8 9 ...
>  $ Subject : Factor w/ 18 levels "308","309","310",..: 1 1 1 1 1 1 1 1 1 1 ...
>  - attr(*, "terms")=Classes 'terms', 'formula' length 3 Reaction ~ Days
>   .. ..- attr(*, "variables")= language list(Reaction, Days)
>   .. ..- attr(*, "factors")= int [1:2, 1] 0 1
>   .. .. ..- attr(*, "dimnames")=List of 2
>   .. .. .. ..$ : chr [1:2] "Reaction" "Days"
>   .. .. .. ..$ : chr "Days"
>   .. ..- attr(*, "term.labels")= chr "Days"
>   .. ..- attr(*, "order")= int 1
>   .. ..- attr(*, "intercept")= int 1
>   .. ..- attr(*, "response")= int 1
>   .. ..- attr(*, ".Environment")=<R_GlobalEnv>
>   .. ..- attr(*, "predvars")= language list(Reaction, Days)
>   .. ..- attr(*, "dataClasses")= Named chr [1:2] "numeric" "numeric"
>   .. .. ..- attr(*, "names")= chr [1:2] "Reaction" "Days"
>  - attr(*, "na.action")=Class 'omit'  Named int [1:3] 14 34 92
>   .. ..- attr(*, "names")= chr [1:3] "14" "34" "92"
> 
> At this point one can get the fitted values with the NA's in the
> correct positions via the thoroughly unobvious
> 
>> stats:::naresid.exclude(attr(fm1 at frame, "na.action"), fm1 at mu)
>   [1] 253.4849 273.1808 292.8766 312.5725 332.2684 351.9642 371.6601 391.3559
>   [9] 411.0518 430.7476 213.1921 214.8561 216.5201       NA 219.8481 221.5121
>  [17] 223.1761 224.8401 226.5041 228.1681 212.8984 217.8564 222.8143 227.7723
>  [25] 232.7302 237.6882 242.6462 247.6041 252.5621 257.5200 275.9163 281.5164
>  [33] 287.1166       NA 298.3168 303.9170 309.5171 315.1173 320.7174 326.3175
>  [41] 273.5320 280.9433 288.3545 295.7658 303.1771 310.5884 317.9997 325.4110
>  [49] 332.8223 340.2336 260.3778 270.5809 280.7841 290.9872 301.1903 311.3934
>  [57] 321.5965 331.7996 342.0027 352.2058 268.1086 278.3689 288.6293 298.8896
>  [65] 309.1500 319.4103 329.6707 339.9310 350.1914 360.4517 244.2257 255.7619
>  [73] 267.2981 278.8342 290.3704 301.9066 313.4427 324.9789 336.5151 348.0512
>  [81] 251.2783 250.9591 250.6398 250.3205 250.0013 249.6820 249.3627 249.0435
>  [89] 248.7242 248.4050 283.8511       NA 322.6857 342.1031 361.5204 380.9377
>  [97] 400.3550 419.7723 439.1897 458.6070 226.4061 238.0217 249.6373 261.2529
> [105] 272.8685 284.4841 296.0997 307.7153 319.3310 330.9466 238.3395 255.4259
> [113] 272.5123 289.5986 306.6850 323.7714 340.8578 357.9442 375.0305 392.1169
> [121] 256.0056 263.4524 270.8991 278.3458 285.7925 293.2393 300.6860 308.1327
> [129] 315.5795 323.0262 272.0277 286.0636 300.0995 314.1354 328.1714 342.2073
> [137] 356.2432 370.2791 384.3150 398.3510 254.6441 265.9888 277.3335 288.6782
> [145] 300.0229 311.3676 322.7123 334.0570 345.4017 356.7464 225.9405 241.2159
> [153] 256.4914 271.7668 287.0422 302.3176 317.5930 332.8684 348.1438 363.4192
> [161] 252.2315 261.7074 271.1833 280.6592 290.1351 299.6110 309.0870 318.5629
> [169] 328.0388 337.5147 263.5955 275.3631 287.1308 298.8985 310.6662 322.4338
> [177] 334.2015 345.9692 357.7368 369.5045
> 
> There is some subtlety about the difference between na.omit and
> na.exclude that eludes me.  I'm sure it is documented and I could,
> with time, unravel the mysteries but I really would prefer not to
> spend that time.  I was hoping that someone else knew what the
> distinction was and how one would construct an extractor that would
> return the desired value and not be such an abominable kludge as the
> above.

Yes. I've prepared and tested the following patch:

Index: C:/Users/hennil/Desktop/R_packages/source/lme4_0.999375-8/R/lmer.R
===================================================================
--- C:/Users/hennil/Desktop/R_packages/source/lme4_0.999375-8/R/lmer.R 
(revision 99)
+++ C:/Users/hennil/Desktop/R_packages/source/lme4_0.999375-8/R/lmer.R 
(working copy)
@@ -882,7 +882,14 @@
        })

  setMethod("fitted", signature(object = "mer"),
-	  function(object, ...) object at eta)
+	  function(object, ...)
+      {
+          eta <- object at eta
+          na.action <- attr(object at frame, "na.action")
+          if(!is.null(na.action))
+              napredict(na.action, eta)
+          else eta
+      })

  setMethod("formula", signature(x = "mer"),
  	  function(x, ...)


After applying the patch, the distinction between `na.omit' and 
`na.exclude' is easily shown:

 > fm1 <- lmer(Reaction ~ Days + (Days|Subject), slps)
 > fitted(fm1)
   [1] 253.4849 273.1808 292.8766 312.5725 332.2684 351.9642
   [7] 371.6601 391.3559 411.0518 430.7476 213.1921 214.8561
  [13] 216.5201 219.8481 221.5121 223.1761 224.8401 226.5041
  [19] 228.1681 212.8984 217.8564 222.8143 227.7723 232.7302
  [25] 237.6882 242.6462 247.6041 252.5621 257.5200 275.9163
  [31] 281.5164 287.1166 298.3168 303.9170 309.5171 315.1173
  [37] 320.7174 326.3175 273.5320 280.9433 288.3545 295.7658
  [43] 303.1771 310.5884 317.9997 325.4110 332.8223 340.2336
  [49] 260.3778 270.5809 280.7841 290.9872 301.1903 311.3934
  [55] 321.5965 331.7996 342.0027 352.2058 268.1086 278.3689
  [61] 288.6293 298.8896 309.1500 319.4103 329.6707 339.9310
  [67] 350.1914 360.4517 244.2257 255.7619 267.2981 278.8342
  [73] 290.3704 301.9066 313.4427 324.9789 336.5151 348.0512
  [79] 251.2783 250.9591 250.6398 250.3205 250.0013 249.6820
  [85] 249.3627 249.0435 248.7242 248.4050 283.8511 322.6857
  [91] 342.1031 361.5204 380.9377 400.3550 419.7723 439.1897
  [97] 458.6070 226.4061 238.0217 249.6373 261.2529 272.8685
[103] 284.4841 296.0997 307.7153 319.3310 330.9466 238.3395
[109] 255.4259 272.5123 289.5986 306.6850 323.7714 340.8578
[115] 357.9442 375.0305 392.1169 256.0056 263.4524 270.8991
[121] 278.3458 285.7925 293.2393 300.6860 308.1327 315.5795
[127] 323.0262 272.0277 286.0636 300.0995 314.1354 328.1714
[133] 342.2073 356.2432 370.2791 384.3150 398.3510 254.6441
[139] 265.9888 277.3335 288.6782 300.0229 311.3676 322.7123
[145] 334.0570 345.4017 356.7464 225.9405 241.2159 256.4914
[151] 271.7668 287.0422 302.3176 317.5930 332.8684 348.1438
[157] 363.4192 252.2315 261.7074 271.1833 280.6592 290.1351
[163] 299.6110 309.0870 318.5629 328.0388 337.5147 263.5955
[169] 275.3631 287.1308 298.8985 310.6662 322.4338 334.2015
[175] 345.9692 357.7368 369.5045
 >
 > fm2 <- lmer(Reaction ~ Days + (Days|Subject), slps,
+             na.action = "na.exclude")
 > fitted(fm2)
   [1] 253.4849 273.1808 292.8766 312.5725 332.2684 351.9642
   [7] 371.6601 391.3559 411.0518 430.7476 213.1921 214.8561
  [13] 216.5201       NA 219.8481 221.5121 223.1761 224.8401
  [19] 226.5041 228.1681 212.8984 217.8564 222.8143 227.7723
  [25] 232.7302 237.6882 242.6462 247.6041 252.5621 257.5200
  [31] 275.9163 281.5164 287.1166       NA 298.3168 303.9170
  [37] 309.5171 315.1173 320.7174 326.3175 273.5320 280.9433
  [43] 288.3545 295.7658 303.1771 310.5884 317.9997 325.4110
  [49] 332.8223 340.2336 260.3778 270.5809 280.7841 290.9872
  [55] 301.1903 311.3934 321.5965 331.7996 342.0027 352.2058
  [61] 268.1086 278.3689 288.6293 298.8896 309.1500 319.4103
  [67] 329.6707 339.9310 350.1914 360.4517 244.2257 255.7619
  [73] 267.2981 278.8342 290.3704 301.9066 313.4427 324.9789
  [79] 336.5151 348.0512 251.2783 250.9591 250.6398 250.3205
  [85] 250.0013 249.6820 249.3627 249.0435 248.7242 248.4050
  [91] 283.8511       NA 322.6857 342.1031 361.5204 380.9377
  [97] 400.3550 419.7723 439.1897 458.6070 226.4061 238.0217
[103] 249.6373 261.2529 272.8685 284.4841 296.0997 307.7153
[109] 319.3310 330.9466 238.3395 255.4259 272.5123 289.5986
[115] 306.6850 323.7714 340.8578 357.9442 375.0305 392.1169
[121] 256.0056 263.4524 270.8991 278.3458 285.7925 293.2393
[127] 300.6860 308.1327 315.5795 323.0262 272.0277 286.0636
[133] 300.0995 314.1354 328.1714 342.2073 356.2432 370.2791
[139] 384.3150 398.3510 254.6441 265.9888 277.3335 288.6782
[145] 300.0229 311.3676 322.7123 334.0570 345.4017 356.7464
[151] 225.9405 241.2159 256.4914 271.7668 287.0422 302.3176
[157] 317.5930 332.8684 348.1438 363.4192 252.2315 261.7074
[163] 271.1833 280.6592 290.1351 299.6110 309.0870 318.5629
[169] 328.0388 337.5147 263.5955 275.3631 287.1308 298.8985
[175] 310.6662 322.4338 334.2015 345.9692 357.7368 369.5045
 >

I hope you find this useful.


Henric



> 
> 
> On Fri, Mar 14, 2008 at 9:43 AM, Douglas Bates <bates at stat.wisc.edu> wrote:
>> On Fri, Mar 14, 2008 at 7:31 AM, Andy Fugard <a.fugard at ed.ac.uk> wrote:
>>
>>  >  Is there any way to get fitted values out of lmer so that any NAs in the
>>  >  predictors/outcome pop out in the correct position?
>>
>>  Yes.
>>
>>
>>  >   From what I can tell, the "fitted" function doesn't have access to
>>  >  missingness in the inputs.
>>
>>  It does not, in its current implementation.  I believe there is
>>  information regarding the pattern of missingness in the model frame
>>  (stored in the "frame" slot of the fitted model) and that can be used
>>  to generate an appropriate fitted response.
>>
>>  Does anyone know exactly how this is done for linear models?  If so,
>>  are you willing to provide a patch to the fitted method for mer
>>  objects?  I can put it on my ToDo list but when the first item on the
>>  list is "finish the book" you almost never get to the second and
>>  subsequent items.
>>
> 
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>




More information about the R-sig-mixed-models mailing list