## [1] "Y" "X1" "X2" "X3" "X4" "X5"
Objetivo: Resolver un problema práctico usando lo aprendido del análisis de regresión lineal. Problema. En un estudio a gran escala realizado en EE.UU sobre la eficacia en el control de infecciones hospitalarias se recogió información en 113 hospitales, los datos se encuentran en publicados junto con el este archivo (datos2.txt). La base de datos contiene las siguientes columnas (variables):
\(Y:\) Riesgo de Infección
\(X1:\) Duración de la estadía
\(X2:\) Rutina de cultivos
\(X3:\) Número de camas
\(X4:\) Censo promedio diario
\(X5:\) Número de enfermeras
Preguntas a resolver.
1. Estime un modelo de regresión lineal múltiple que explique el Riesgo de Infección en términos de todas las variables predictoras. Analice la significancia de la regresión y de los parámetros individuales. Interprete los parámetros estimados. Calcule e interprete el coeficiente de determinación múltiple R2. Comente los resultados.
##
## Call:
## lm(formula = Y ~ ., data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.03944 -0.80043 -0.00266 0.60450 2.23292
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.5986009 1.5159559 -0.395 0.694365
## X1 0.2106683 0.0785765 2.681 0.009501 **
## X2 0.0197512 0.0277108 0.713 0.478803
## X3 0.0470925 0.0132888 3.544 0.000779 ***
## X4 0.0105604 0.0073166 1.443 0.154213
## X5 0.0008996 0.0007379 1.219 0.227679
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.035 on 59 degrees of freedom
## Multiple R-squared: 0.4468, Adjusted R-squared: 0.3999
## F-statistic: 9.53 on 5 and 59 DF, p-value: 1.058e-06
De acuerdo al resultado anterior el modelo ajustado que explica el Riesgo de InfeccIón en términos de todas las variables predictoras es el siguiente: \[{Y}_i={\beta_0}+{\beta_1}X_{i1}+{\beta_2}X_{i2}+{\beta_3}X_{i3}+{\beta_4}X_{i4}+{\beta_5}X_{i5}+\epsilon_i, \; i=1,2,...,65.\]
\[{Y}_i=-0.5986+0.2107X_{i1}+0.0198X_{i2}+0.0471X_{i3}+0.0106X_{i4}+0.0009X_{i5}+ \epsilon_i , \; i=1,2,...,65.\]
Podemos ver en el summary del modelo u \(p-value: 1.058e-06\) y si manejamos un nivel de significancia del 0.05 (\(\alpha = 0.05\)), podemos concluir que el modelo es significativo y al menos una de las variables es valida para explicar el modelo, por lo tanto es útil para explicar el riesgo de infección.
Prueba de hipótesis:
\(H_0:\beta_j=0\) vs \(H_1:\beta_j\neq0\) con \(j=1,..,5\) con un \(\alpha=0.05\)
Criterio de Rechazo:
\(H_0\) si el \(valor-p= p(t_{\alpha/2,n-p}>|t_0|)\) es pequeño
Con un \(\alpha<0.05\) y \(\beta_1 = 0.009501\ _y \ \beta_3=0.000779\) podemos concluir que son significativos ya que el \(\alpha > \beta_1, \beta_3\), lo cual nos indica que se rechaza la hipótesis nula. A la luz del problema lo que nos quiere decir es que la Duración de la estadía y Número de camas nos ayuda a explicar el Riesgo de Infección en los hospitales.
Las otras variables no cumplen el supuesto por ende podemos concluir que no son significativas.
Se tiene un \(R^2=0.4468\) lo que nos indica que el \(44.68%\) de la variabilidad total es explicada por el modelo. En otras palabras variabilidad presente en el riesgo de infrección.
2. Use la tabla de todas las regresiones posibles, para probar la significancia simultánea del subconjunto de tres variables con los valores p mayores del punto anterior. Según el resultado de la prueba es posible descartar del modelo las variables del subconjunto?.
##
## Call:
## lm(formula = Y ~ X2 + X4 + X5, data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.9012 -0.8825 -0.1692 0.5834 3.0467
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.1711795 1.6920689 0.692 0.491464
## X2 0.0110235 0.0296634 0.372 0.711463
## X4 0.0269514 0.0074407 3.622 0.000596 ***
## X5 0.0016713 0.0008099 2.063 0.043329 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.196 on 61 degrees of freedom
## Multiple R-squared: 0.2357, Adjusted R-squared: 0.1981
## F-statistic: 6.27 on 3 and 61 DF, p-value: 0.0008843
El resultado de la prueba nos indica que no es posible descartar todas las variables del subconjunto, ya que las variables X5 y X4, bajo esté modelo ajustado, son significativas para explicar la variable de respuesta.
3. Plantee una pregunta donde su solución implique el uso exclusivo de una prueba de hipótesis lineal general de la forma \(H_0: L\beta = 0\) (solo se puede usar este procedimiento y no SSextra), donde especifique claramente la matriz L, el modelo reducido y la expresión para el estadístico de prueba.
\[{Y}_i=\beta_0+\beta_1X_{i1}+\beta_2X_{i2}+\beta_3X_{i3}+\beta_4X_{i4}+\beta_5X_{i5}+ \epsilon_i, \; \epsilon_i\thicksim^{iid}N(0,\sigma^2)\;\; i=1,2,...,65.\]
Vamos a probar la siguiente hipótesis:
Igualamos las ecuaciones a cero y tenemos:
Estamos probando \(m=3\) ecuaciones diferentes, entonce armamos una matriz \(L\) de orden \(m*p\) en este caso \(3*6\)
\[L = \begin{equation} \begin{pmatrix} \beta_0 & \beta_1 & \beta_2 & \beta_3 & \beta_4 & \beta_5 \\ 0 & 1 & -1 & 0 & 0 & 0\\ 0 & 0 & 0 & 1 & 0 & -1\\ 0 & 0 & 0 & 0 & 1 & 0 \end{pmatrix} \end{equation} \begin{pmatrix} \beta_0\\ \beta_1\\ \beta_2\\ \beta_3\\ \beta_4\\ \beta_5 \end{pmatrix} = \begin{pmatrix} 0\\ 0\\ 0\\ \end{pmatrix}\]
Se plantea el modelo reducido bajo la hipótesis nula
Modelo Reducido \[{Y}_i=\beta_0+\beta_1X_{i1}+\beta_1X_{i2}+\beta_3X_{i3}+\beta_3X_{i5}+\epsilon_i\]
\[Y_i= \beta_0+\beta_1(X_{i1}+X_{i2})+\beta_3(X_{i3}+X_{i5})+\epsilon_i\]
## Linear hypothesis test
##
## Hypothesis:
## X1 - X2 = 0
## X3 - X5 = 0
## X4 = 0
##
## Model 1: restricted model
## Model 2: Y ~ X1 + X2 + X3 + X4 + X5
##
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 62 100.530
## 2 59 63.203 3 37.327 11.615 4.399e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Con un \(alpha= 0.05\) y un \(p-value= 4.399e-06\) podemos concluir que \(\alpha>p-value\), por lo tanto hay evidencia suficiente para rechazar hipótesis nula, lo que nos indica que el Censo promedio diario es significativa y nos ayuda a explicar el Riesgo de infección en los hospitales dado que las demás variables están en el modelo.
4. Realice una validación de los supuestos en los errores y examine si hay valores atípicos, de balanceo e influenciables. Qué puede decir acerca de la validez de éste modelo?. Argumente.
## Test stat Pr(>|Test stat|)
## X1 -0.4564 0.64980
## X2 -0.1161 0.90797
## X3 -1.5808 0.11935
## X4 0.2852 0.77647
## X5 -2.1379 0.03675 *
## Tukey test -1.8235 0.06823 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Shapiro-Wilk normality test
##
## data: restud
## W = 0.98663, p-value = 0.7094
Primero hay que tener en cuenta que los outliers son los datos se encuentran alejados del eje \(Y\) y los de balanceo son los que se encuentran alejados del eje \(X\).
En la gráfica de \(residuales \ vs \ X1\) se puede apreciar un dato que posiblemente sea de balanceo ya que se encuentra muy alejado de los demás a la derecha, posiblemente influenciable, también podemos apreciar que en la gráfica \(residuales \ vs \ X3\) hay un dato en -1.5 en \(Y\) y 61 en \(X\), que también puede ser de balanceo e influenciable.
En el gráfico de residuales vs Normal, vemos que los datos se encuentran ajustados a la recta normal, con cola a la izquierda pesada y posiblemente un dato atípico.
influence.measures(modelo)
## Influence measures of
## lm(formula = Y ~ ., data = datos) :
##
## dfb.1_ dfb.X1 dfb.X2 dfb.X3 dfb.X4 dfb.X5 dffit cov.r
## 1 2.99e-02 7.91e-02 -0.06191 -0.041093 -0.027867 1.81e-02 -0.116785 1.176
## 2 1.41e-01 -3.15e-02 -0.19008 -0.085804 0.089953 1.56e-01 -0.336245 0.920
## 3 1.47e-02 -1.79e-02 -0.00972 -0.023011 0.020129 1.38e-02 0.045178 1.142
## 4 6.03e-02 2.32e-02 0.00935 0.065691 -0.296731 3.83e-01 0.552310 0.883
## 5 4.23e-02 9.56e-02 -0.02736 0.056667 -0.216995 1.97e-01 0.344268 1.143
## 6 1.08e-02 -2.22e-02 -0.00890 -0.050711 0.031524 6.35e-02 0.093682 1.171
## 7 3.04e-02 -1.83e-02 -0.17885 0.070859 0.264872 1.48e-01 -0.485371 0.918
## 8 1.01e-02 1.64e-03 -0.01123 -0.010044 0.001261 2.25e-03 -0.015533 1.200
## 9 1.49e-02 2.81e-02 -0.01328 -0.030262 -0.018352 -1.52e-02 0.056633 1.151
## 10 3.15e-02 7.56e-03 -0.03932 -0.047879 0.032579 -1.06e-03 0.069235 1.159
## 11 3.19e-03 -1.49e-02 -0.00632 -0.008700 0.038508 6.52e-03 0.064820 1.116
## 12 1.07e-01 -6.30e-02 -0.06206 0.005880 -0.035033 4.27e-02 0.148990 1.134
## 13 7.20e-03 -1.56e-01 0.08863 -0.024046 0.048639 -1.01e-01 0.313609 0.894
## 14 -3.07e-02 -1.89e-02 0.03697 0.001626 0.002423 3.96e-02 0.053984 1.303
## 15 3.52e-02 -5.96e-02 -0.00442 0.005699 -0.027241 2.27e-02 -0.115185 1.099
## 16 4.94e-03 1.69e-01 -0.00281 0.049096 -0.196029 -7.40e-02 0.251440 1.150
## 17 -7.22e-02 2.50e-02 0.05855 0.087741 -0.070884 7.55e-02 -0.175259 1.078
## 18 9.57e-03 1.16e-03 -0.01259 -0.009385 0.010937 -3.43e-03 0.017907 1.257
## 19 -1.27e-01 1.24e-01 -0.00682 -0.029917 0.129416 7.92e-02 -0.304631 1.027
## 20 7.52e-02 2.40e-01 -0.15867 0.002108 -0.149775 3.12e-02 -0.329959 1.093
## 21 -4.36e-02 -1.45e-01 0.02273 -0.044362 0.295269 -3.84e-02 0.336039 1.134
## 22 4.41e-02 1.38e-01 -0.07501 -0.000197 -0.018245 -3.54e-01 -0.370128 1.310
## 23 -2.29e-01 5.11e-02 0.15432 -0.038838 0.155118 2.81e-02 0.289280 1.386
## 24 1.63e-02 -4.62e-03 -0.02164 0.009829 0.036085 -3.34e-02 0.081163 1.124
## 25 2.60e-02 6.20e-03 -0.03406 -0.215578 0.047595 5.41e-02 -0.250729 1.112
## 26 3.11e-02 -9.51e-03 0.00360 0.007094 -0.100104 1.11e-01 0.165476 1.280
## 27 2.07e-01 -5.17e-01 0.01299 0.085322 0.017165 1.30e-01 -0.571053 1.923
## 28 -2.38e-01 -1.52e-01 0.25029 0.182674 0.118516 -7.76e-02 -0.367513 1.238
## 29 1.98e-01 -1.99e-01 -0.10512 -0.054937 0.094779 -8.63e-02 0.338128 1.031
## 30 1.85e-01 1.90e-01 -0.09872 -0.073106 -0.321658 -5.15e-01 -0.658574 1.164
## 31 -6.96e-02 5.08e-02 0.06250 0.042064 -0.038010 -1.29e-02 0.122865 1.117
## 32 -2.43e-01 -2.18e-01 0.30004 0.196432 0.014296 1.55e-01 -0.381692 1.164
## 33 -3.91e-02 4.84e-03 -0.01551 -0.119847 0.150939 1.26e-01 0.235368 1.237
## 34 -1.06e-04 -1.29e-04 0.00025 0.000105 -0.000274 -3.45e-05 -0.000536 1.155
## 35 -1.10e-01 7.65e-03 0.09842 0.068765 -0.061433 1.23e-01 -0.208335 1.071
## 36 3.68e-02 -2.83e-02 -0.02157 0.011111 -0.024702 2.30e-02 -0.076008 1.139
## 37 1.60e-01 2.29e-01 -0.29164 -0.377534 0.045854 1.51e-01 -0.567576 0.766
## 38 2.12e-02 -4.72e-02 -0.03892 0.009793 0.072427 4.42e-02 -0.121168 1.165
## 39 3.21e-02 1.23e-02 -0.00456 0.226807 -0.122886 -8.68e-02 0.268598 1.253
## 40 -1.12e-02 2.38e-02 0.04659 -0.128583 -0.106477 5.02e-02 -0.251068 1.136
## 41 -3.98e-01 -3.86e-01 0.52220 -0.078449 0.322114 -6.74e-02 0.742942 0.974
## 42 8.64e-03 2.39e-02 0.04456 -0.088463 -0.057398 -3.57e-02 0.322918 0.683
## 43 -1.46e-01 -1.44e-01 0.23861 0.345880 -0.065172 -1.44e-01 0.473715 0.925
## 44 -1.09e-01 7.84e-05 0.11141 0.045808 -0.024929 1.02e-02 -0.122716 1.315
## 45 1.45e-01 2.39e-01 -0.28280 0.511462 -0.047121 -2.53e-01 0.784904 1.181
## 46 -6.30e-02 1.17e-01 0.00538 0.079436 0.038603 -1.36e-01 0.258468 1.029
## 47 -2.09e-02 -5.89e-03 0.02904 0.011557 -0.006889 -5.93e-03 0.042201 1.141
## 48 -1.55e-02 -2.46e-03 -0.00285 0.005891 0.032577 2.49e-02 -0.057703 1.165
## 49 -1.04e-02 7.17e-03 0.00423 -0.010188 0.003268 6.19e-03 -0.030041 1.129
## 50 -2.50e-02 1.98e-01 -0.04522 0.041391 -0.055807 -2.23e-01 -0.299688 1.072
## 51 1.02e-01 4.40e-02 0.03152 -0.067397 -0.293871 -9.03e-02 0.422717 0.942
## 52 -1.78e-03 2.65e-02 0.03662 -1.154133 0.392677 -4.99e-01 -1.402146 1.157
## 53 4.34e-05 3.44e-03 -0.00184 -0.003332 0.001948 -3.35e-04 0.006031 1.171
## 54 -9.33e-02 2.72e-02 0.02911 0.061865 0.078624 2.37e-02 -0.180324 1.111
## 55 -2.70e-03 2.86e-02 -0.07367 0.118066 0.026206 1.26e-01 -0.315231 0.912
## 56 6.45e-02 1.10e-01 -0.06422 0.029968 -0.147023 4.40e-02 0.219403 1.090
## 57 -2.02e-01 9.00e-02 0.03938 0.083172 0.206564 6.92e-02 -0.399519 1.010
## 58 7.06e-03 -1.13e-02 -0.00236 -0.012475 -0.005845 2.39e-02 -0.040781 1.147
## 59 1.60e-02 -3.63e-02 -0.01445 0.028631 0.021665 -8.05e-04 -0.073947 1.141
## 60 1.88e-01 -2.07e-01 0.00936 -0.001215 -0.173935 5.26e-03 0.453769 0.737
## 61 -1.42e-02 1.72e-02 0.01014 0.015831 -0.023994 6.40e-03 -0.040532 1.153
## 62 -1.87e-02 4.89e-03 0.06605 0.040664 -0.124709 -1.31e-01 -0.250528 1.021
## 63 4.66e-02 -2.00e-02 -0.04440 -0.002789 0.027755 -1.76e-02 0.066652 1.275
## 64 1.53e-01 -9.27e-02 -0.19628 0.052581 0.146516 6.75e-02 -0.334150 1.179
## 65 4.51e-02 1.81e-02 -0.03924 -0.032042 -0.022137 7.71e-02 0.177881 0.988
## cook.d hat inf
## 1 2.31e-03 0.0740
## 2 1.84e-02 0.0468
## 3 3.46e-04 0.0350
## 4 4.90e-02 0.0879
## 5 1.98e-02 0.1161
## 6 1.48e-03 0.0660
## 7 3.82e-02 0.0804
## 8 4.09e-05 0.0772
## 9 5.43e-04 0.0445
## 10 8.11e-04 0.0525
## 11 7.10e-04 0.0245
## 12 3.74e-03 0.0583
## 13 1.60e-02 0.0378
## 14 4.94e-04 0.1510
## 15 2.23e-03 0.0322
## 16 1.06e-02 0.0947
## 17 5.15e-03 0.0427
## 18 5.44e-05 0.1187
## 19 1.54e-02 0.0630
## 20 1.81e-02 0.0921
## 21 1.88e-02 0.1101
## 22 2.30e-02 0.2007 *
## 23 1.41e-02 0.2242 *
## 24 1.11e-03 0.0335
## 25 1.05e-02 0.0772
## 26 4.63e-03 0.1486
## 27 5.49e-02 0.4471 *
## 28 2.26e-02 0.1655
## 29 1.89e-02 0.0733
## 30 7.14e-02 0.2014
## 31 2.54e-03 0.0425
## 32 2.43e-02 0.1352
## 33 9.34e-03 0.1362
## 34 4.87e-08 0.0407
## 35 7.25e-03 0.0497
## 36 9.77e-04 0.0405
## 37 5.08e-02 0.0681
## 38 2.48e-03 0.0686
## 39 1.21e-02 0.1516
## 40 1.06e-02 0.0882
## 41 8.90e-02 0.1567
## 42 1.63e-02 0.0201 *
## 43 3.64e-02 0.0793
## 44 2.55e-03 0.1639 *
## 45 1.01e-01 0.2349
## 46 1.11e-02 0.0513
## 47 3.02e-04 0.0343
## 48 5.64e-04 0.0546
## 49 1.53e-04 0.0228
## 50 1.49e-02 0.0759
## 51 2.91e-02 0.0713
## 52 3.13e-01 0.3432 *
## 53 6.17e-06 0.0540
## 54 5.46e-03 0.0567
## 55 1.62e-02 0.0408
## 56 8.06e-03 0.0595
## 57 2.63e-02 0.0832
## 58 2.82e-04 0.0379
## 59 9.25e-04 0.0416
## 60 3.24e-02 0.0435
## 61 2.78e-04 0.0426
## 62 1.04e-02 0.0470
## 63 7.53e-04 0.1335
## 64 1.87e-02 0.1298
## 65 5.24e-03 0.0230
Verificando en las otras medidas como cov.r podemos ver que los datos que son influenciables son el 22, 23, 27, 44, 52.
Para verificar las observaciones de balanceo tiene que cumplir la formula: \(\frac{2(k+1)}{n}= \frac{2(5+1)}{65}= 0.184\), verificando en la tabla los “hat inf”, los números superiores a \(0.184\). Tenemos que observaciones de balanceo son: 22,23,27,30,45.
5. Verificar la presencia de multicolinealidad usando graficos y/o indicadores apropiados.
vif(modelo)
## X1 X2 X3 X4 X5
## 1.375666 1.176970 1.270949 1.302309 1.124675
colldiag <- function(mod,scale=TRUE,center=FALSE,add.intercept=TRUE) {
result <- NULL
if (center) add.intercept<-FALSE
if (is.matrix(mod)||is.data.frame(mod)) {
X<-as.matrix(mod)
nms<-colnames(mod)
}
else if (!is.null(mod$call$formula)) {
X<-mod$model[,-1] # delete the dependent variable
}
X<-na.omit(X) # delete missing cases
if (add.intercept) {
X<-cbind(1,X) # add the intercept
colnames(X)[1]<-"intercept"
}
X<-scale(X,scale=scale,center=center)
svdX<-svd(X)
svdX$d
condindx<-svdX$d[1]/svdX$d
Phi=svdX$v%*%diag(1/svdX$d)
Phi<-t(Phi^2)
pi<-prop.table(Phi,2)
dim(condindx)<-c(length(condindx),1)
colnames(condindx)<-"cond.index"
rownames(condindx)<-1:nrow(condindx)
colnames(pi)<-colnames(X)
result$condindx<-condindx
result$pi<-pi
class(result)<-"colldiag"
result
}
colldiag(modelo)
## $condindx
## cond.index
## 1 1.000000
## 2 4.225132
## 3 4.778789
## 4 12.542973
## 5 16.120536
## 6 37.127569
##
## $pi
## intercept X1 X2 X3 X4 X5
## [1,] 0.0002386136 0.0008820943 0.0002496554 0.006424605 0.001403204 0.007901015
## [2,] 0.0001855077 0.0001147921 0.0001862061 0.155728657 0.003353482 0.776605732
## [3,] 0.0023903752 0.0049986896 0.0038354112 0.640772446 0.004296292 0.124324670
## [4,] 0.0223457381 0.0037191274 0.0285460562 0.114363084 0.920850823 0.002073082
## [5,] 0.0399258551 0.9761466034 0.0293340090 0.004189488 0.065410293 0.067694494
## [6,] 0.9349139102 0.0141386933 0.9378486621 0.078521720 0.004685906 0.021401005
##
## attr(,"class")
## [1] "colldiag"
Para el análisis que haremos a continuación debemos tener muy presente que si el \(vif > 10\) tenemos problema de multicolinealidad. Como ninguna variable nos muestra que superior a 10 entonces no se detecta problema de multicolinealidad, ahora evaluemos la proporción de varianzas.
R calcula la raíz cuadrada del numero de condición \(\sqrt{k_6}= 37.128\) lo que detecta multicolinealidad debido a que supero a 31. El resto de las columnas son la proporción de varianzas y como bien sabemos si supera el 0.5 nos dice que hay una relación lineal entre dos covariables. Se mira por filas. Como no se aprecia ningún VIF que supere el 0.5 no se puede concluir que no hay problemas de multicolinealidad.