Regresión lineal con restricción de pendiente

Question

Quiero realizar una regresión lineal muy sencilla en R . La fórmula es tan sencilla como $y = ax + b$ . Sin embargo, me gustaría que la pendiente ( $a$ ) esté dentro de un intervalo, digamos, entre 1,4 y 1,6.

¿Cómo se puede hacer esto?

Answer 1

Quiero realizar una regresión lineal en R. ... Me gustaría que la pendiente estuviera dentro de un intervalo, digamos, entre 1,4 y 1,6. ¿Cómo se puede hacer esto?

(i) De forma sencilla:

• encajar la regresión. Si está dentro de los límites, has terminado.

• Si no está dentro de los límites, establece la pendiente al límite más cercano, y

• estimar el intercepto como la media de $(y - ax)$ sobre todas las observaciones.

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

Edición: en realidad (como se menciona en el ejemplo de abajo), en vanilla R, nls puede hacer restricciones de caja; como se muestra en el ejemplo, eso es realmente muy fácil de hacer.

Se puede utilizar la regresión restringida más directamente; creo que el pcls del paquete "mgcv" y la función nnls del paquete "nnls".

--

Editar para responder a la pregunta de seguimiento -

Iba a mostrarte cómo usarlo con nlminb ya que eso viene con R, pero me di cuenta de que nls ya utiliza las mismas rutinas (las rutinas PORT) para implementar mínimos cuadrados restringidos, así que mi ejemplo de abajo hace ese caso.

Nota: en mi ejemplo de abajo, $a$ es el intercepto y $b$ es la pendiente (la convención más común en estadística). Me he dado cuenta después de ponerlo aquí de que has empezado al revés; no obstante, voy a dejar el ejemplo "al revés" en relación con tu pregunta.

En primer lugar, 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 de ella, sólo causada por la variación aleatoria. Así que vamos a utilizar 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 ves, tienes una pendiente justo en el límite. Si pasas el modelo ajustado a summary incluso producirá errores estándar y valores t, pero no estoy seguro de que sean significativos/interpretables.

Entonces, ¿cómo se compara mi sugerencia (1)? (es decir, ajustar la pendiente al límite más cercano y promediar los residuos $y-bx$ para estimar el intercepto)

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

Es la misma estimación ...

En el gráfico siguiente, la línea azul son los mínimos cuadrados y la línea roja son los mínimos cuadrados restringidos:

Answer 2

El segundo método de Glen_b, que utiliza los mínimos cuadrados con una restricción de caja, puede aplicarse más fácilmente mediante la regresión de crestas. La solución de 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). Así que, siguiendo la sugerencia de Whuber, 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.

La ventaja de este enfoque es que no se necesitan herramientas de optimización sofisticadas, sólo la regresión de crestas, que ya está disponible en R (y en muchos otros paquetes).

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

Answer 3

Otro enfoque sería utilizar métodos bayesianos para ajustar la regresión y elegir una distribución a priori sobre $a$ que sólo tiene soporte en la región que desea, por ejemplo, un uniforme de 1,4 a 1,6, o una distribución beta desplazada y escalada a ese dominio.

Hay muchos ejemplos en la web y en los programas informáticos sobre el uso de métodos bayesianos para la regresión, basta con seguir uno de esos ejemplos y cambiar la prioridad en $a$ .

Este resultado seguirá dando 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).

Answer 4

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

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