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, para este caso analizaremos una muestra aleatoria de 65 hospitales, donde contaremos con las siguientes variables

Abreviatura y nombre de la Variable Descripción
Y: Riesgo de infección Probabilidad promedio estimada de adquirir infección en el hospital (en porcentaje).
X1: Duración de la estadía Duración promedio de la estadía de todos los pacientes en el hospital (en días).
X2: Rutina de cultivos Razón del número de cultivos realizados en pacientes sin síntomas de infección hospitalaria, por cada 100.
X3: Número de camas Número promedio de camas en el hospital durante el periodo del estudio.
X4: Censo promedio diario Número promedio de pacientes en el hospital por día durante el periodo del estudio.
X5:Número de enfermeras Número promedio de enfermeras, equivalentes a tiempo completo, durante el periodo del estudio.

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

library(readr)
Hospitales <- read_csv("Hospitales.csv")
## Rows: 65 Columns: 6
## -- Column specification --------------------------------------------------------
## Delimiter: ","
## dbl (6): Y, X1, X2, X3, X4, X5
## 
## i Use `spec()` to retrieve the full column specification for this data.
## i Specify the column types or set `show_col_types = FALSE` to quiet this message.
View(Hospitales)
attach(Hospitales)
myPanel.cor=function(x, y, digits = 2, prefix = "", cex.cor){
  usr <- par("usr"); on.exit(par(usr = usr))
  par(usr = c(0, 1, 0, 1))
  r <- cor(x, y)
  txt <- format(c(r, 0.123456789), digits = digits)[1]
  txt <- paste(prefix, txt, sep = "")
  if(missing(cex.cor))
    cex = 0.4/strwidth(txt)
  text(0.5, 0.5, txt, cex = 1 + 1.5*abs(r))
}
myPanel.box=function(x, ...){
  usr <- par("usr", bty = 'n')
  on.exit(par(usr))
  par(usr = c(-1, 1, min(x) - 0.5, max(x) + 0.5))
  b <- boxplot(x, plot = F)
  whisker.i <- b$stats[1,]
  whisker.s <- b$stats[5,]
  hinge.i <- b$stats[2,]
  mediana <- b$stats[3,]
  hinge.s <- b$stats[4,]
  rect(-0.5, hinge.i, 0.5, mediana, col = 'gray')
  segments(0, hinge.i, 0, whisker.i, lty = 2)
  segments(-0.1, whisker.i, 0.1, whisker.i)
  rect(-0.5, mediana, 0.5, hinge.s, col = 'gray')
  segments(0, hinge.s, 0, whisker.s, lty = 2)
  segments(-0.1, whisker.s, 0.1, whisker.s)
}
pairs(Hospitales, lower.panel = myPanel.cor,upper.panel = panel.smooth, diag.panel = myPanel.box, labels = names(Hospitales))

Según el análisis de la matriz de gráficos de dispersión anterior podemos definir nuestro modelo de la siguiente manera

\[Y_i = \beta_0 + \beta_1x_{i1} + \beta_2x_{i2} + \beta_3x_{i3} + \beta_4x_{i4} + \beta_5x_{i5} + \varepsilon_i,\quad i= 1, 2, \ldots, 65\]

con los siguientes supuestos: \[\varepsilon_i \overset{\text{iid}}{\sim} N\left(0,\sigma ^2 \right), \quad i=1, 2, \ldots, 65\]

#Definición modelo de regresión LM, Tabala ANOVA y Tabla de parámetros estimados.
modelo <- lm(Y ~ ., Hospitales)
anova(modelo)
## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value    Pr(>F)    
## X1         1 27.872 27.8721 26.0185 3.770e-06 ***
## X2         1  0.499  0.4995  0.4662    0.4974    
## X3         1 19.028 19.0275 17.7621 8.703e-05 ***
## X4         1  2.054  2.0542  1.9176    0.1713    
## X5         1  1.592  1.5919  1.4861    0.2277    
## Residuals 59 63.203  1.0712                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(modelo)$coefficients
##                  Estimate   Std. Error    t value     Pr(>|t|)
## (Intercept) -0.5986008691 1.5159559354 -0.3948669 0.6943650156
## X1           0.2106683091 0.0785764898  2.6810603 0.0095010402
## X2           0.0197512078 0.0277108065  0.7127619 0.4788026393
## X3           0.0470924905 0.0132888409  3.5437621 0.0007788529
## X4           0.0105603543 0.0073165927  1.4433432 0.1542126323
## X5           0.0008995893 0.0007379464  1.2190442 0.2276789260

Teniendo estos valores podemos construir nuestro modelo estimado: \[\widehat{Y}=-0.598600869+0.2106683091X_{i1}+0.0197512078X_{i2}+0.0470924905X_{i3}+0.0105603543X_{i4}+ 0.0008995893X_{i5}\]

Significancia de los parámetros Se establece el siguiente juego de hipótesis:

\[\begin{array}{1} H_0: \beta_j = 0\\ H_1: \beta_j \ne 0 \end{array}\ \text { con }\ j = 0, 1,\ldots, 5.\] A través de los resultados del estadístico t y el valor p, obtenidos en las dos columnas de la tabla de coeficientes estimados, a un nivel de significancia de \(\alpha\) = \(0.05\) se puede concluir que los parámetros individuales \(\widehat{\beta}_1\) \(\widehat{\beta}_3\) son significativos en presencia de los demás parámetros, por otro lado, los otros parámetros no significativos.

Interpretación de los parámetros \(\widehat{\beta}_1\): Indica que por cada unidad que aumente el número de días en el hospital (\(X_3\)), el promedio de riesgo de infección disminuye en promedio 0,2106683091 unidades cuando las demás predictoras se mantienen constantes \(\widehat{\beta}_3\): Indica que por cada unidad que aumente el número de camas (\(X_3\)) el promedio del riesgo de infeccion disminuye en promedio 0.0470924905 unidades cuando las demas predictoras se mantengas constantes.

Significancia de la regresión se establece el siguiente juego de hipótesis:

\[\begin{aligned} H_0: &\ \beta_1 = \beta_2 = \cdots = \beta_5 = 0 \quad \text { vs } \\H_1: &\ \text { algún } \beta_j \quad \neq 0, j = 1, \ldots, 5 \end{aligned}\]

Para lograr determinar la significancia de la regresión se utiliza la tabla de análisis de varianza ANOVA

myAnova=function(modeloreg){
SSq<-unlist(anova(modeloreg)["Sum Sq"])
k<-length(SSq)-1
SSR<-sum(SSq[1:k])
SSE<-SSq[(k+1)]
MSR<-SSR/k
df.error<-unlist(anova(modeloreg)["Df"])[k+1]
MSE<-SSE/df.error
F0<-MSR/MSE
VP<-pf(F0,k,df.error,lower.tail=F)
result<-
data.frame(SumSq=c(SSR,SSE),Df=c(k,df.error),MeanSq=c(MSR
,MSE),
F0=c(round(F0,digits=3),' 
'),P.value=c(format(VP,scientific = TRUE,digits=3),' '),
row.names =c("Modelo","Error"))
cat("Tabla ANOVA Modelo de Regresi�n","\n")
result
}

myAnova(modelo)
## Tabla ANOVA Modelo de Regresi<U+FFFD>n
##           SumSq Df    MeanSq   F0  P.value
## Modelo 51.04527  5 10.209053 9.53 1.06e-06
## Error  63.20335 59  1.071243   \n

De la tabla ANOVA anterior se obtienen los valores del estadistico de prueba F = 9.53 y su correspondiente valor p = 1.06e-06 Como el valor p es menor a 0.05 = \(\alpha\) se rechaza la hipotesis nula, es decir, el modelo de regresion lineal multiple es significativo. Por lo que el control del riesgo de in feccion depende significativamente de alguna de las varibles del modelo.

Interpretacion del \(R^2\)

El \(R^2\) se obtiene de la división entre el SSR y el SST, de la tabla ANOVA encontramos que el SSR = 51.04527 y el SST = SSR + SSE, por tanto, el SST = 114.24862. Al hacer la división tenemos que \(R^2\) = 0.4467 y en porcentaje 44.67% El 44.67% de la variabilidad total en el control del riesgo de contagio es explicado por el modelo de regresión lineal múltiple propuesto

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?.

Para ello debemos usar la tabla de todas las regresiones posibles ya que, para poder probar la significancia de las variables, debemos probar la significancia simultánea de sus respectivos parámetros, equivaliéndose a la siguiente prueba de hipótesis:

\[\begin{aligned} H_0: &\ \beta_1 = \beta_3 \text { vs } \\H_1: &\ \text { algún } \beta_j \neq 0, j = 1,3. \end{aligned}\] Para probar esta hipótesis se usan las sumas de cuadrados extra, para ello definimos los modelos: Completo y reducido (bajo 𝐻0) Para conocer los SSR Y SSE de cada modelo. Para ello usamos la funcion MyAllRegTable() la cual nos traerá una tabla con todas las posibles regresiones y variables con intercepto disponibles. La tabla nos da {\(2^K\)} − 1 Número de modelos, siendo 𝑘 el número de variables predictoras, en este caso 5, lo que nos trae 31 modelos.

myAllRegTable=function(lm.model, response = model.response(model.frame(lm.model)), MSE = F){
  regTable <- summary(regsubsets(model.matrix(lm.model)[, -1], response,
                                 nbest = 2^(lm.model$rank - 1) - 1, really.big = T))
  pvCount <- as.vector(apply(regTable$which[, -1], 1, sum))
  pvIDs <- apply(regTable$which[, -1], 1, function(x) as.character(paste(colnames(model.matrix(lm.model)[, -1])[x],
                                                                         collapse = " ")))
  result <- if(MSE){
    data.frame(k = pvCount, R_sq = round(regTable$rsq, 3), adj_R_sq = round(regTable$adjr2, 3),
               MSE = round(regTable$rss/(nrow(model.matrix(lm.model)[,-1]) - (pvCount + 1)), 3),
               Cp = round(regTable$cp, 3), Variables_in_model = pvIDs)
  } else {
    data.frame(k = pvCount, R_sq = round(regTable$rsq, 3), adj_R_sq = round(regTable$adjr2, 3),
               SSE = round(regTable$rss, 3),
               Cp = round(regTable$cp, 3), Variables_in_model = pvIDs)
  }
  format(result, digits = 6)
}
library(leaps)
myAllRegTable(modelo)
##    k  R_sq adj_R_sq     SSE     Cp Variables_in_model
## 1  1 0.249    0.237  85.813 19.106                 X3
## 2  1 0.244    0.232  86.376 19.632                 X1
## 3  1 0.182    0.169  93.468 26.252                 X4
## 4  1 0.068    0.053 106.508 38.425                 X5
## 5  1 0.001   -0.014 114.090 45.503                 X2
## 6  2 0.411    0.392  67.257  3.784              X1 X3
## 7  2 0.319    0.297  77.804 13.629              X3 X4
## 8  2 0.311    0.289  78.739 14.503              X1 X4
## 9  2 0.293    0.270  80.764 16.392              X3 X5
## 10 2 0.276    0.253  82.667 18.169              X2 X3
## 11 2 0.258    0.235  84.721 20.087              X1 X5
## 12 2 0.248    0.224  85.877 21.166              X1 X2
## 13 2 0.234    0.209  87.519 22.698              X4 X5
## 14 2 0.182    0.156  93.416 28.203              X2 X4
## 15 2 0.071    0.041 106.102 40.046              X2 X5
## 16 3 0.430    0.402  65.108  3.778           X1 X3 X4
## 17 3 0.422    0.393  66.088  4.693           X1 X3 X5
## 18 3 0.415    0.386  66.850  5.404           X1 X2 X3
## 19 3 0.359    0.327  73.279 11.405           X3 X4 X5
## 20 3 0.336    0.303  75.868 13.822           X2 X3 X4
## 21 3 0.328    0.295  76.791 14.684           X1 X4 X5
## 22 3 0.325    0.292  77.090 14.963           X2 X3 X5
## 23 3 0.314    0.280  78.400 16.186           X1 X2 X4
## 24 3 0.261    0.224  84.459 21.842           X1 X2 X5
## 25 3 0.236    0.198  87.321 24.514           X2 X4 X5
## 26 4 0.442    0.405  63.748  4.508        X1 X3 X4 X5
## 27 4 0.433    0.395  64.795  5.486        X1 X2 X3 X4
## 28 4 0.427    0.389  65.435  6.083        X1 X2 X3 X5
## 29 4 0.379    0.338  70.904 11.188        X2 X3 X4 X5
## 30 4 0.329    0.284  76.656 16.558        X1 X2 X4 X5
## 31 5 0.447    0.400  63.203  6.000     X1 X2 X3 X4 X5

Para la prueba de Hipótesis tenemos cómo estadístico de prueba a: \[F_0 = \dfrac{\text{MSextra}}{\text{MSE}} = MSR(\beta_2,\beta_4,\beta_5|\beta_0,\beta_1,\beta_3)/{\text{MSE}}= \dfrac{\text{(SSE}{(\beta_0,\beta_1,\beta_3}) - \text{SSE}(\beta_0,\beta_2,\beta_4,\beta_5)/)r}{\text{MSE}}\] \[((67.257-63.203)/3)/1.071 =1.26\] Para el criterio de decisión necesitamos el valor crítico de una distribución \(f_{3,59}\)a un nivel de significancia de 𝛼 = 0,05 Esto es \(f_{0.05,3,59}\) =2.760767

Como nuestro estadistico es menor que la f, No podemos rechazar \(H_0\) y concluimos que es posible descartar estas variables del modelo, aunque los valores sean muy cercanos, ya que el riesgo de contagio no depende de estas variables

3. Plantee una pregunta donde su solución implique el uso exclusivo de una prueba de hipótesis lineal general de la forma H0 : L = 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.

¿Es la rutina de cultivos, el censo promedio diario, y el número promedo de enfermeras significantes a igual nivel? ¿Es el número de camas en el hospital es igual de significante que el número de días promedio en el hospital?. Considere la prueba de significancia del modelo de RLM con K= 5 variables predictoras y el siguiente juego de hipótesis:

\[\begin{aligned} H_0: &\ \beta_2 = \beta_4 = \beta_5, \quad \beta_1 = \beta_3 \quad \text { vs }\\ H_1: &\ \ \beta_2 \neq \beta_4 \neq \beta_4, \quad \beta_1 \neq \beta_3 \end{aligned}\]

Podemos reescribir la hipotesis nula de la siguiente manera:

\[\begin{aligned} H_0: &\ \beta_2 - \beta_4=0,\quad \beta_2 - \beta_4 = 0,\quad \beta_2 - \beta_5=0, \quad \beta_1 - \beta_3 = 0 \end{aligned}\] De manera que la hipótesis nula contiene m = 4 ecuaciones y se puede escribir como:

\[H_0: \left\{\begin{array}{c} \beta_2 - \beta_4 = 0\\ \beta_2 - \beta_5 = 0\\ \beta_4 - \beta_5 = 0\\ \beta_1 - \beta_3 = 0 \end{array}\right.\]

De forma matricial:

\[H_0: \begin{array}{cccccc111} \left[\begin{array}{ccccccccc} 0 & 0 & 1 & 0 & -1 & 0\\ 0 & 0 & 1 & 0 & 0 & -1\\ 0 & 0 & 0 & 0 & 1 & -1\\ 0 & 1 & 0 & -1 & 0 & 0 \end{array}\right] & \left[\begin{array}{c} \beta_0\\ \beta_1\\ \beta_2\\ \beta_3\\ \beta_4\\ \beta_5 \end{array}\right] & = & \left[\begin{array}{c} 0\\ 0\\ 0\\ 0 \end{array}\right] \end{array}\]

Por lo tanto la hipótesis lineal general con r = 4 filas linealmente independientes es:

\[\boldsymbol{L}: \begin{array}{cccccc111} \left[\begin{array}{ccccccccc} 0 & 0 & 1 & 0 & -1 & 0\\ 0 & 0 & 1 & 0 & 0 & -1\\ 0 & 0 & 0 & 0 & 1 & -1\\ 0 & 1 & 0 & -1 & 0 & 0 \end{array}\right] & \left[\begin{array}{c} \beta_0\\ \beta_1\\ \beta_2\\ \beta_3\\ \beta_4\\ \beta_5 \end{array}\right] & = & \left[\begin{array}{c} 0\\ 0\\ 0\\ 0 \end{array}\right] \end{array}\]

A partir de esto nuestro modelo reducido será

\[\begin{aligned} \boldsymbol{RM}: \boldsymbol{Y} = \beta_0 + \beta_1(X_1+X_3) + \beta_2(X_2+X_4+X_5)+\beta_3(X_4+X_5) + \varepsilon\\ \\\end{aligned}\] Definimos: \[F_0 = \dfrac{\text{MSH}}{\text{MSE}} = \dfrac{\text{SSH/}r}{\text{MSE}}= \dfrac{[\text{SSE}(\text{RM}) = \text{SSE}(\text{FM})]/r}{\text{MSE}}\] 4. Realice una validación de los supuestos en los errores y examine si hay valores atípicos, de balanceo einfluenciales. Qué puede decir acerca de la validez de éste modelo?. Argumente.

Validación supuestos de normalidad y varianza constante

#Gráfica y prueba normalidad Shapiro-wilk
myQQnorm=function(modelo, student = F, ...){
  if(student){
    res <- rstandard(modelo)
    lab.plot <- "Normal Q-Q Plot of Studentized Residuals"
  } else {
    res <- residuals(modelo)
    lab.plot <- "Normal Q-Q Plot of Residuals"
  }
  shapiro <- shapiro.test(res)
  shapvalue <- ifelse(shapiro$p.value < 0.001, "P value < 0.001", paste("P value = ", round(shapiro$p.value, 4), sep = ""))
  shapstat <- paste("W = ", round(shapiro$statistic, 4), sep = "")
  q <- qqnorm(res, plot.it = FALSE)
  qqnorm(res, main = lab.plot, ...)
  qqline(res, lty = 2, col = 2)
  text(min(q$x, na.rm = TRUE), max(q$y, na.rm = TRUE)*0.95, pos = 4, 'Shapiro-Wilk Test', col = "blue", font = 2)
  text(min(q$x, na.rm = TRUE), max(q$y, na.rm = TRUE)*0.80, pos = 4, shapstat, col = "blue", font = 3)
  text(min(q$x, na.rm = TRUE), max(q$y, na.rm = TRUE)*0.65, pos = 4, shapvalue, col = "blue", font = 3)
}

myQQnorm(modelo)

La recta punteada representa un ajuste de la distribución normal, el patrón de los residuales sigue esta recta, se concluye por el motivo anterior y apoyado por la prueba de normalidad Shapiro - Wilks indica que los errores son normales, el supuesto de normalidad se cumple

Supuesto de varianza constante - Grafica de residuales vs valores ajustados

res.stud= round(rstandard(modelo))
yhat=round(modelo$fitted.values)
plot(yhat, res.stud, xlab = "Valores Ajustados", ylab = "Residuales Estudentizados")
abline(h = 0, lty = 2, lwd = 2, col = 2)

Deseamos Probar:

\[H_0: \ V[\varepsilon_i] = \sigma^2 \ \text{ vs. } \ \ H_1: \ V[\varepsilon_i]\neq \sigma^2\] Al hacer el análisis de la gráfica se nota que los puntos estan distribuidos uniformemente o centrados al inicio ,sin embargo se nota una disperción importante al final, se concluye que el supuesto de varianza constante no se cumple.

## Diagnósticos para identificar valores extremos

# Cálculo de errores estándar de los valores ajustados
se.yhat <- round(predict(modelo, se.fit = T)$se.fit, 4)
# Residuales crudos del modelo
residuals <- round(modelo$residuals, 4)
# Distancias de Cook
Cooks.D <- round(cooks.distance(modelo), 4)
# Valores de la diagonal de la matriz H
hii.value <- round(hatvalues(modelo), 4)
# Dffits
Dffits <- round(dffits(modelo), 4)
# Tabla de diagnósticos
data.frame(Hospitales, yhat, se.yhat, residuals, res.stud, Cooks.D, hii.value, Dffits)
##      Y    X1   X2   X3    X4  X5 yhat se.yhat residuals res.stud Cooks.D
## 1  3.7  7.58 56.7 20.8  88.0  97    4  0.2816   -0.4143        0  0.0023
## 2  2.8  9.97 58.2 16.5  76.5  90    4  0.2239   -1.5171       -2  0.0184
## 3  4.2  8.88 51.5 10.1  86.9 305    4  0.1938    0.2430        0  0.0003
## 4  6.2 10.15 51.9 16.4  59.2 568    4  0.3069    1.7268        2  0.0490
## 5  5.7 11.18 51.0 18.8  55.9 595    5  0.3527    0.9251        1  0.0198
## 6  4.5  9.61 52.4  6.9  87.2 487    4  0.2659    0.3552        0  0.0015
## 7  1.6  8.82 58.2  3.8  51.7  80    3  0.2935   -1.6059       -2  0.0382
## 8  5.1 10.30 59.6 27.8  88.9 175    5  0.2875   -0.0539        0  0.0000
## 9  4.1 10.47 53.2  5.7  69.1 196    4  0.2183    0.2677        0  0.0005
## 10 4.4 10.02 49.5  8.3  93.0 265    4  0.2371    0.2986        0  0.0008
## 11 5.0  9.78 52.3 17.6  95.9 270    5  0.1621    0.4208        0  0.0007
## 12 4.3  7.65 47.1 16.4  65.7 318    4  0.2500    0.6045        1  0.0037
## 13 5.3  8.15 54.9 12.3  79.8  99    4  0.2012    1.5863        2  0.0160
## 14 4.8  9.84 62.2 12.0  82.3 600    5  0.4022    0.1231        0  0.0005
## 15 4.4 11.65 54.5 18.6  96.1 248    5  0.1858   -0.6460       -1  0.0022
## 16 5.3 11.77 54.1 17.3  56.0 196    5  0.3185    0.7681        1  0.0106
## 17 2.9  8.86 51.3  9.5  87.5 100    4  0.2139   -0.8425       -1  0.0051
## 18 4.3  9.89 45.2 11.8 108.7 190    4  0.3566    0.0478        0  0.0001
## 19 2.0  7.08 52.0 12.3  56.4  87    3  0.2598   -1.1731       -1  0.0154
## 20 2.7  7.14 57.6 13.1  92.6  92    4  0.3142   -1.0208       -1  0.0181
## 21 5.6  8.95 53.7 18.9 122.8 147    5  0.3435    0.9334        1  0.0188
## 22 4.1  9.35 53.8 15.9  80.9 833    5  0.4636   -0.6862       -1  0.0230
## 23 6.6 13.95 65.9 15.6 133.5 356    6  0.4901    0.4935        1  0.0141
## 24 5.1  9.76 50.9 21.9  97.0 150    5  0.1895    0.4465        0  0.0011
## 25 4.5 10.05 52.0 36.7  87.5 184    5  0.2876   -0.8635       -1  0.0105
## 26 4.3  9.23 51.6 11.6  42.6 620    4  0.3989    0.3811        0  0.0046
## 27 6.5 19.56 59.9 17.2 113.7 306    7  0.6921   -0.4911       -1  0.0549
## 28 2.9 10.79 44.2  2.6  56.6 461    4  0.4211   -0.7824       -1  0.0226
## 29 4.5  6.70 48.6 13.0  80.8  76    3  0.2802    1.1934        1  0.0189
## 30 4.9 11.07 53.2 28.5 122.0 768    6  0.4645   -1.2056       -1  0.0714
## 31 5.6 11.48 57.6 20.3  82.0 252    5  0.2134    0.5938        1  0.0025
## 32 3.0 11.20 45.0  7.0  78.9 130    4  0.3806   -0.9295       -1  0.0243
## 33 5.7 11.80 53.8  9.1 116.9 571    5  0.3820    0.5734        1  0.0093
## 34 5.0 11.03 49.9 19.7 102.1 318    5  0.2088   -0.0027        0  0.0000
## 35 2.9  8.90 49.7 12.7  86.9  52    4  0.2307   -0.9205       -1  0.0073
## 36 4.5 11.46 56.9 15.6  97.7 191    5  0.2084   -0.3777        0  0.0010
## 37 2.5  8.54 56.1 27.0  82.5  98    5  0.2701   -2.0394       -2  0.0508
## 38 3.4 10.42 58.0  8.0  59.0 119    4  0.2711   -0.4490        0  0.0025
## 39 5.8  9.50 49.3 42.0  70.9  98    5  0.4030    0.6087        1  0.0121
## 40 4.8 10.24 49.0 36.3 112.6 195    6  0.3073   -0.8004       -1  0.0106
## 41 5.4  7.93 64.1  7.5  98.1  68    4  0.4097    1.6116        2  0.0890
## 42 6.3  9.74 54.4 11.4  76.1 221    4  0.1468    2.2329        2  0.0163
## 43 6.3  8.84 56.3 29.6  82.6  85    5  0.2915    1.5816        2  0.0364
## 44 3.4  8.45 38.8 12.9  85.0 235    4  0.4190   -0.2644        0  0.0025
## 45 7.8 12.07 43.7 52.4 105.3 157    7  0.5016    1.2718        1  0.1010
## 46 6.4 11.62 53.9 25.5  99.2 133    5  0.2345    1.1180        1  0.0111
## 47 4.6  9.68 57.8 16.7  79.0 186    4  0.1916    0.2297        0  0.0003
## 48 3.1  8.63 54.0  8.4  56.2  76    3  0.2420   -0.2435        0  0.0006
## 49 4.1  9.05 51.2 20.5  79.8 195    4  0.1564   -0.2027        0  0.0002
## 50 2.9  7.91 52.8 11.9  79.5 477    4  0.2851   -1.0397       -1  0.0149
## 51 4.7  8.77 54.5  5.2  47.0 143    3  0.2764    1.5047        2  0.0291
## 52 5.4 11.18 45.7 60.5  85.8 640    7  0.6063   -1.5902       -2  0.3130
## 53 4.8 12.01 52.8 10.8  96.9 298    5  0.2405    0.0256        0  0.0000
## 54 2.3  7.95 51.8  4.6  54.9 163    3  0.2464   -0.7423       -1  0.0055
## 55 2.0  8.93 56.0  6.2  72.5  95    4  0.2091   -1.5318       -2  0.0162
## 56 5.5 11.08 50.2 18.6  63.6 387    5  0.2525    0.8772        1  0.0081
## 57 1.4  7.14 51.7  4.1  45.7 115    3  0.2986   -1.3058       -1  0.0263
## 58 4.7 10.72 53.8 23.2  94.1 113    5  0.2015   -0.2103        0  0.0003
## 59 3.9 11.15 56.5  7.7  73.9 281    4  0.2112   -0.3621        0  0.0009
## 60 5.5  7.63 52.1 11.6  61.1 197    3  0.2157    2.0934        2  0.0324
## 61 3.7  8.48 51.1 12.1  92.8 166    4  0.2135   -0.1963        0  0.0003
## 62 3.9 10.73 50.6 19.3 101.0 445    5  0.2244   -1.1371       -1  0.0104
## 63 4.2  7.53 42.0 23.1  98.9  95    4  0.3781    0.1650        0  0.0008
## 64 2.9 10.80 63.9  1.6  57.4 130    4  0.3729   -0.8372       -1  0.0187
## 65 5.6 10.12 51.7 14.9  79.1 362    4  0.1569    1.1828        1  0.0052
##    hii.value  Dffits
## 1     0.0740 -0.1168
## 2     0.0468 -0.3362
## 3     0.0350  0.0452
## 4     0.0879  0.5523
## 5     0.1161  0.3443
## 6     0.0660  0.0937
## 7     0.0804 -0.4854
## 8     0.0772 -0.0155
## 9     0.0445  0.0566
## 10    0.0525  0.0692
## 11    0.0245  0.0648
## 12    0.0583  0.1490
## 13    0.0378  0.3136
## 14    0.1510  0.0540
## 15    0.0322 -0.1152
## 16    0.0947  0.2514
## 17    0.0427 -0.1753
## 18    0.1187  0.0179
## 19    0.0630 -0.3046
## 20    0.0921 -0.3300
## 21    0.1101  0.3360
## 22    0.2007 -0.3701
## 23    0.2242  0.2893
## 24    0.0335  0.0812
## 25    0.0772 -0.2507
## 26    0.1486  0.1655
## 27    0.4471 -0.5711
## 28    0.1655 -0.3675
## 29    0.0733  0.3381
## 30    0.2014 -0.6586
## 31    0.0425  0.1229
## 32    0.1352 -0.3817
## 33    0.1362  0.2354
## 34    0.0407 -0.0005
## 35    0.0497 -0.2083
## 36    0.0405 -0.0760
## 37    0.0681 -0.5676
## 38    0.0686 -0.1212
## 39    0.1516  0.2686
## 40    0.0882 -0.2511
## 41    0.1567  0.7429
## 42    0.0201  0.3229
## 43    0.0793  0.4737
## 44    0.1639 -0.1227
## 45    0.2349  0.7849
## 46    0.0513  0.2585
## 47    0.0343  0.0422
## 48    0.0546 -0.0577
## 49    0.0228 -0.0300
## 50    0.0759 -0.2997
## 51    0.0713  0.4227
## 52    0.3432 -1.4021
## 53    0.0540  0.0060
## 54    0.0567 -0.1803
## 55    0.0408 -0.3152
## 56    0.0595  0.2194
## 57    0.0832 -0.3995
## 58    0.0379 -0.0408
## 59    0.0416 -0.0739
## 60    0.0435  0.4538
## 61    0.0426 -0.0405
## 62    0.0470 -0.2505
## 63    0.1335  0.0667
## 64    0.1298 -0.3342
## 65    0.0230  0.1779

1. Identificación de observaciones atípicas. Criterio: Se considera que una observación es atípica cuando su residual estudentizado |r_i| es mayor que 3, |r_i|> 3

Según la columna res.stud de la tabla anterior, de valores atípicos, se concluye que no hay valores atípicos

¨2. Identificación de puntos de balanceo. Criterio: La observación 𝑖 es un punto de balanceo si \(h_{ii}\) > 2𝑝/𝑛, donde 𝑝 son los parámetros y 𝑛 el numero de observaciones En este ejercicio se tiene que: \(ℎ_{ii}\) > 2(6/65) = 0,184615 Segun la columna hii.value se obtiene que las observaciones:22,23,24,27,30,45,52. son puntos de balanceo

3. Identificación de observaciones influenciales. Criterios: -La observación 𝑖 sera influencible si \(D_i\) > 1 -Una observación será influencial si \(|DFFITS_i|\) > 2*\(Raiz(p/n)\) En este ejercicio deben superar el valor 0,6076 De acuerdo con la columna Cooks.D de distancias de Cook tenemos que no se muestran influenciales. Sin embargo de acuerdo con la columna Dffits de la tabla de valores DFFITS tenemos que las observaciones 30,41,45, y 52 son influenciales

Conclusiones: • No hay valores atipicos • Las observaciones 22,23,24,27,30, 45 y 52 son puntos de balanceo. • Las observaciones 30,41,45 y 52 son influenciales. Los puntos de balanceo y los influenciales pueden estar afectando el análisis de normalidad de los errores y el supuesto de varianza constante, podemos imaginar que en ausencia de estos el comportamiento de los residuales sería mejor, sin embargo, esto no es motivo para retirarlos del modelo.

5. Verificar la presencia de multicolinealidad usando graficos y/o indicadores apropiados.

cor(Hospitales)
##            Y        X1          X2         X3         X4          X5
## Y  1.0000000 0.4939233  0.03723050  0.4988930 0.42648615  0.26029349
## X1 0.4939233 1.0000000  0.20636081  0.1982719 0.37908182  0.29406162
## X2 0.0372305 0.2063608  1.00000000 -0.2476300 0.03743146 -0.08511824
## X3 0.4988930 0.1982719 -0.24762997  1.0000000 0.35967155  0.10258439
## X4 0.4264861 0.3790818  0.03743146  0.3596716 1.00000000  0.07684465
## X5 0.2602935 0.2940616 -0.08511824  0.1025844 0.07684465  1.00000000

Al no tener, ningún valor mayor a 0.5, hasta el momento no tengo indicios de multicoleanidad, sin embargo la matrix de correlación de variables, no es lo suficientemente fuerte para descartar la multicoleanidad.

library(car)
## Loading required package: carData
vifs=c(0,vif(modelo))
miscoeficientes=function(modeloreg,datosreg)
{coefi=coef(modeloreg)
datos2=as.data.frame(scale(datosreg))
coef.std=c(0,coef(lm(update(formula(modeloreg),~.+0),datos2)))
limites=confint(modeloreg,level=0.95)
vifs=c(0,vif(modeloreg))
resul=data.frame(Estimación=coefi,Limites=limites,Vif=vifs,Coef.Std=coef.std)
cat("Coeficientes estimados, sus I.C, Vifs y Coeficientes estimados 
estandarizados","\n")
resul}
miscoeficientes(modelo,Hospitales)
## Coeficientes estimados, sus I.C, Vifs y Coeficientes estimados 
## estandarizados
##                Estimación Limites.2.5.. Limites.97.5..      Vif   Coef.Std
## (Intercept) -0.5986008691  -3.632021689    2.434819951 0.000000 0.00000000
## X1           0.2106683091   0.053437116    0.367899502 1.375666 0.30449605
## X2           0.0197512078  -0.035697988    0.075200403 1.176970 0.07487653
## X3           0.0470924905   0.020501581    0.073683400 1.270949 0.38685413
## X4           0.0105603543  -0.004080114    0.025200822 1.302309 0.15949429
## X5           0.0008995893  -0.000577038    0.002376217 1.124675 0.12518477

Podemos concluir que bajo el criterio de los factores de inflación de varianza no se han encontrado problemas de multicolinealidad; pero sigue siendo necesario evaluar por demás pruebas para descartar la misma

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
}

misDiagnostcolin=function(modeloreg,centrar=
F){
if(centrar==F){
X=model.matrix(modeloreg)
val.prop=prcomp(X,center=FALSE,scale=TRUE)$sdev^2
Ind=colldiag(modeloreg)
resul=data.frame(Val.propio=val.prop,Ind.Cond=Ind$condindx,Pi=Ind$pi)
cat("Diagnósticos Multicolinealidad -
Intercepto incluído","\n",
"Índices de Condición y Proporciones de 
Varianza","\n")
}
else{
X=model.matrix(modeloreg)[,-1]
val.prop=prcomp(X,center=TRUE,scale=TRUE)$sdev^2
Ind=colldiag(modeloreg,center=TRUE,scale=TRUE)
resul=data.frame(Val.propio=val.prop,Ind.Cond=Ind$condindx,Pi=Ind$pi)
cat("Diagnósticos Multicolinealidad -
Intercepto ajustado","\n",
"Índices de Condición y Proporciones de 
Varianza","\n")
}
resul}
misDiagnostcolin(modelo)
## Diagnósticos Multicolinealidad -
## Intercepto incluído 
##  Índices de Condición y Proporciones de 
## Varianza
##   Val.propio cond.index Pi.intercept        Pi.X1        Pi.X2       Pi.X3
## 1 5.40182563   1.000000 0.0002386136 0.0008820943 0.0002496554 0.006424605
## 2 0.30259372   4.225132 0.0001855077 0.0001147921 0.0001862061 0.155728657
## 3 0.23654019   4.778789 0.0023903752 0.0049986896 0.0038354112 0.640772446
## 4 0.03433520  12.542973 0.0223457381 0.0037191274 0.0285460562 0.114363084
## 5 0.02078651  16.120536 0.0399258551 0.9761466034 0.0293340090 0.004189488
## 6 0.00391875  37.127569 0.9349139102 0.0141386933 0.9378486621 0.078521720
##         Pi.X4       Pi.X5
## 1 0.001403204 0.007901015
## 2 0.003353482 0.776605732
## 3 0.004296292 0.124324670
## 4 0.920850823 0.002073082
## 5 0.065410293 0.067694494
## 6 0.004685906 0.021401005

Al analizar los indices de condición y proporción de varianza con intecerto, podemos inferir la presencia de multicoleanidad moderada para los valores propios 4 y 5. Y una multicoleanidad Fuerte para el Valor propio 6; sin embargo al analizar los \(π_{ij}}\)>0.5, no encontramos multicoleanidad entre dos variables Analicemos los indices de condición y proporción de varianza sin intecerto

misDiagnostcolin(modelo, centrar = T )
## Diagnósticos Multicolinealidad -
## Intercepto ajustado 
##  Índices de Condición y Proporciones de 
## Varianza
##   Val.propio cond.index       Pi.X1        Pi.X2      Pi.X3       Pi.X4
## 1  1.7323654   1.000000 0.132673811 0.0001653996 0.10741655 0.142120663
## 2  1.2396841   1.182127 0.084395758 0.4433688520 0.13165390 0.001072787
## 3  0.9874882   1.324506 0.009577422 0.0282670390 0.05662206 0.136118366
## 4  0.5461281   1.781035 0.013604334 0.2065052224 0.64568619 0.572581500
## 5  0.4943342   1.872015 0.759748675 0.3216934870 0.05862130 0.148106685
##          Pi.X5
## 1 6.484536e-02
## 2 2.244896e-06
## 3 6.375120e-01
## 4 2.662746e-03
## 5 2.949776e-01

Al analizar los indices de condición y proporción de varianza sin intecerto, podemos inferir la presencia de multicoleanidad moderada para el valor propi 4.al analizar los \(π_{ij}}\)>0.5, no encontramos multicoleanidad entre dos variables, vemos que las variables \(X_3\) y \(X_4\), estan generando una multicoleanidad moderada