El argumento de los pesos en lm y lme es muy diferente en R- ¿Los estoy usando correctamente?

8

Entonces, me parece que la función de pesos en lm da más peso a las observaciones cuanto mayor sea el valor de "peso" de la observación asociada, mientras que la función lme en lme hace exactamente lo contrario. Esto se puede verificar con una simple simulación.

#make 3 vectors- c is used as an uninformative random effect for the lme model
a<-c(1:10)
b<-c(2,4,6,8,10,100,14,16,18,20)
c<-c(1,1,1,1,1,1,1,1,1,1)

Si ahora ejecutara un modelo donde ponderara las observaciones basadas en el inverso de la variable dependiente en lm, solo podría generar exactamente el mismo resultado en nlme si pesa solo con la variable dependiente, sin tomar el inverso.

summary(lm(b~a,weights=1/b))
summary(lme(b~a,random=~1|c,weights=~b))

Puede voltear esto y ver que lo contrario es cierto: especificar pesos = b en lm requiere pesos = 1 / b para obtener un resultado lme coincidente.

Entonces, entiendo esto, solo quiero validar una cosa y hacer una pregunta sobre otra.

  1. Si quiero ponderar mis datos basados ​​en el inverso de la variable dependiente, ¿está bien codificar pesos = ~ (variable dependiente) dentro de lme?
  2. ¿Por qué se escribe lme para manejar pesos completamente diferentes a lm? ¿Cuál es el propósito de esto aparte de generar confusión?

¡Cualquier idea sería apreciada!

colin
fuente
1
La respuesta a 2. es que fueron escritos por personas muy diferentes para hacer cosas muy diferentes. lm()necesitaba ser compatible con S y varios libros, nlme no lo hizo, y pretendía ser más flexible, permitiendo que la heterogeneidad se modelara de manera más flexible de lo que lmpermite.
Gavin Simpson

Respuestas:

12

Q1

En lmela notación weights = ~ b, la varFixedfunción de varianza se usaría con un único argumento b. Esta función agregaría al modelo una función de varianza que tiene la forma, donde toma los valores del argumento vector .s2(v)s2(v)=El |vEl |vb

Por lo tanto, se debe utilizar weights = ~ I(1/b)en lme()tener la varianza de .εyo=1/ /siyo

En lmlo que pasa weightsparece ser exactamente lo contrario; weightses inversamente proporcional a la varianza.

No estoy 100% seguro de lo que quiere decir por peso de mis datos , pero si se refiere a proporcionar la variación heterogénea de las observaciones, entonces creo que quiere weights = ~ I(1/b).

Q2

Mi primera impresión (que tendría que preguntar a los respectivos autores de las dos funciones) es que esto es por culpa lm()y lme()fueron escritos por personas muy diferentes de hacer las cosas muy diferentes. lm()necesitaba (se deseaba que fuera) para ser compatible con S y varios libros, nlme no, y pretendía ser más flexible, permitiendo modelar la heterogeneidad de manera más flexible de lo que lmpermite mediante el uso de funciones de varianza a través de la varFuncinfraestructura.

Gavin Simpson
fuente
Esto es lo suficientemente claro. Por 'ponderar mis datos' quiero decir que quiero que el ajuste del modelo tenga en cuenta que se deben esperar grandes residuales de las grandes observaciones, y que se ajuste a un porcentaje de mínimos cuadrados, en lugar de mínimos cuadrados ordinarios. TAMBIÉN: eliminé la publicación cruzada en el desbordamiento de pila, lo siento
colin
Es posible que desee ver otras funciones de varianza en nlme entonces. Lo que está haciendo es decir que las variaciones de sus observaciones son exactamente el valor (absoluto) de b. Parece mejor decir que la varianza aumentó con b. varPower()por ejemplo tendría la varianza comoσ^2×El |siEl |2δ con δestimó un parámetro modelo. Esto está bien si bno toma 0 valores. Si puede tomar valores 0, entonces la varExp()función puede ser mejor, allí la varianza esvunar(εyo)=σ^2×mi2δ×siyo.
Gavin Simpson
En lm(), tenga en cuenta la redacción de que la varianza es proporcional a la inversa de weights. En el lmecódigo que discutimos, b es la varianza. Siguiendo su explicación, no creo que realmente quiera eso ... También tenga en cuenta que si la varianza aumenta con la respuesta media, entonces un GLMM puede ser apropiado y el paquete lme4 sería adecuado, ya que puede modelar la relación media-varianza directamente , en lugar de mediante la modificación de la matriz de covarianza, que es lo lmeque está haciendo el código.
Gavin Simpson
Finalmente, lo siento si sonaba malhumorado en Stack Overflow . No fue intencional. Olvidé que no puedes votar para cerrar como OT y migrar a Cross Validated . Tienes que dejar un comentario sobre por qué, pero ya había dejado el primer comentario. No elija un sitio SE para su pregunta en función de la cantidad de ojos que lo verán. Elige el lugar más apropiado. No hay nada de malo en promocionar su pregunta en Cross Validated para obtener más ojos, incluso puede publicar el enlace en la sala de chat pública R en Stack Overflow . Crossposting o publicar preguntas OT diluye el recurso si tenemos demasiados, por lo tanto, cerca califican etc.
Gavin Simpson