Regresión lineal con restricción de pendiente

18

Quiero realizar una regresión lineal muy simple en R. La fórmula es tan simple como . Sin embargo, me gustaría que la pendiente ( ) esté dentro de un intervalo, digamos, entre 1.4 y 1.6.y=unX+siun

¿Cómo se puede hacer esto?

Iñigo Hernáez Corres
fuente

Respuestas:

24

Quiero realizar ... regresión lineal en R. ... Me gustaría que la pendiente esté dentro de un intervalo, digamos, entre 1.4 y 1.6. ¿Cómo se puede hacer esto?

(i) Manera simple:

  • ajustarse a la regresión. Si está en los límites, ya está.

  • Si no está en los límites, establezca la pendiente en el límite más cercano y

  • estimar la intersección como el promedio de sobre todas las observaciones.(y-unX)

(ii) Forma más compleja: hacer mínimos cuadrados con restricciones de caja en la pendiente; muchas rutinas de optimización implementan restricciones de cuadro, por ejemplo nlminb(que viene con R) lo hace.

Editar: en realidad (como se menciona en el ejemplo a continuación), en vanilla R, nlspuede hacer restricciones de cuadro; como se muestra en el ejemplo, eso es realmente muy fácil de hacer.

Puede usar la regresión restringida más directamente; Creo que la pclsfunción del paquete "mgcv" y la nnlsfunción del paquete "nnls" sí lo hacen.

-

Editar para responder la pregunta de seguimiento -

Te iba a mostrar cómo usarlo nlminbya que viene con R, pero me di cuenta de que nlsya usa las mismas rutinas (las rutinas PORT) para implementar mínimos cuadrados restringidos, por lo que mi ejemplo a continuación hace ese caso.

NB: en mi ejemplo a continuación, es la intersección yb es la pendiente (la convención más común en las estadísticas). Me di cuenta después de ponerlo aquí que comenzaste al revés; Sin embargo, voy a dejar el ejemplo "hacia atrás" en relación con su pregunta.unsi

Primero, configure algunos datos con la pendiente 'verdadera' dentro del rango:

 set.seed(seed=439812L)
 x=runif(35,10,30)
 y = 5.8 + 1.53*x + rnorm(35,s=5)  # population slope is in range
 plot(x,y)
 lm(y~x)

Call:
lm(formula = y ~ x)

Coefficients:
(Intercept)            x  
     12.681        1.217  

... pero la estimación de LS está muy por fuera, solo causada por una variación aleatoria. Entonces usemos la regresión restringida en nls:

 nls(y~a+b*x,algorithm="port",
   start=c(a=0,b=1.5),lower=c(a=-Inf,b=1.4),upper=c(a=Inf,b=1.6))

Nonlinear regression model
  model: y ~ a + b * x
   data: parent.frame()
    a     b 
9.019 1.400 
 residual sum-of-squares: 706.2

Algorithm "port", convergence message: both X-convergence and relative convergence (5)

Como puede ver, obtiene una pendiente justo en el límite. Si le pasa el modelo ajustado summary, incluso producirá errores estándar y valores t, pero no estoy seguro de cuán significativos / interpretables sean.

y-siX

 b=1.4
 c(a=mean(y-x*b),b=b)
       a        b 
9.019376 1.400000

Es la misma estimación ...

En la gráfica a continuación, la línea azul son los mínimos cuadrados y la línea roja son los mínimos cuadrados restringidos:

línea restringida y LS

Glen_b -Reinstate a Monica
fuente
Gracias por esta respuesta, pero ... ¿podría dar un ejemplo usando alguna de estas funciones?
Iñigo Hernáez Corres
1
+1 Encontrar intervalos de confianza en las estimaciones de parámetros será un desafío en cualquier caso.
whuber
@ IñigoHernáezCorres ve la actualización de mi respuesta, donde ilustre usarlo nlspara hacerlo.
Glen_b: reinstala a Monica el
¡+1 gran respuesta con conexiones sobre dos formas de hacerlo!
Haitao Du
15

El segundo método de Glen_b, el uso de mínimos cuadrados con una restricción de cuadro puede implementarse más fácilmente mediante regresión de cresta. La solución a la regresión de cresta puede verse como el lagrangiano para una regresión con un límite en la magnitud de la norma del vector de peso (y, por lo tanto, su pendiente). Entonces, siguiendo la sugerencia de Whuber a continuación, el enfoque sería restar una tendencia de (1.6 + 1.4) / 2 = 1.5 y luego aplicar la regresión de cresta y aumentar gradualmente el parámetro de cresta hasta que la magnitud de la pendiente sea menor o igual a 0.1.

El beneficio de este enfoque es que no se requieren herramientas de optimización sofisticadas, solo regresión de cresta, que ya está disponible en R (y muchos otros paquetes).

Sin embargo, la solución simple de Glen_b (i) me parece sensata (+1)

Dikran Marsupial
fuente
55
Esto es inteligente, pero ¿estás seguro de que funcionará como se describe? Me parece que el enfoque apropiado sería eliminar una tendencia de (1.6 + 1.4) / 2 = 1.5 y luego controlar el parámetro de la cresta hasta que el valor absoluto de la pendiente sea menor o igual a 0.1.
whuber
1
Sí, esa es una mejor sugerencia. El enfoque de regresión de cresta es realmente más apropiado si la restricción está en la magnitud de la pendiente, ¡parece un problema bastante extraño! Mi respuesta se inspiró originalmente en el comentario de Glen_b sobre las restricciones de cuadro, la regresión de cresta es básicamente una forma más fácil de implementar las restricciones de cuadro.
Dikran Marsupial
Aunque aprecio su reconocimiento de mis comentarios, en mi opinión, esto distrae del contenido de su respuesta. Estamos todos juntos en esto para mejorar nuestro trabajo siempre que podamos, por lo que es suficiente reconocimiento de que actuó de acuerdo con mis sugerencias. Por eso te mereces el incremento de reputación. Si se mueve para realizar ediciones adicionales, considere simplificar el texto eliminando ese material superfluo.
whuber
Se edita material superfluo, sin embargo, disfruto de las colaboraciones y siempre trato de darles a los colaboradores el crédito que se merecen, y sigo pensando moralmente que mereces la mitad de los votos positivos. ; o)
Dikran Marsupial
10

un

un

Este resultado aún dará intervalos creíbles de los parámetros de interés (por supuesto, el significado de estos intervalos se basará en la razonabilidad de su información previa sobre la pendiente).

Greg Snow
fuente
+1, este fue mi primer pensamiento también. Me gustan las otras sugerencias, pero esta me parece la mejor.
gung - Restablece a Monica
0

Otro enfoque podría ser reformular su regresión como un problema de optimización y usar un optimizador. No estoy seguro de si se puede reformular de esta manera, pero pensé en esta pregunta cuando leí esta publicación de blog sobre optimizadores R:

http://zoonek.free.fr/blosxom/R/2012-06-01_Optimization.html

Wayne
fuente