MASTER UNIVERSITARIO IDeA

Profesor: Francisco Pla Martos

Departamento de Matemáticas.
Facultad de Ciencias y Tecnologías Químicas. Ciudad Real.
Universidad de Castilla-La Mancha, UCLM

Autor: Francisco Plá Martos 1

 MASTER UNIVERSITARIO IDeA

T7.2 Regresión Lineal Simple y

Múltiple
REGRESIÓN LINEAL SIMPLE

Ejemplo1:

Los datos sobre la longitud de las plantas en centímetros al cabo de un

año de vida son:

x = (15.3, 17.8, 20.7, 25.1, 16.4, 21.6, 19.6, 18.8, 20.2, 19.4)

y los datos de las mismas plantas cuando son adultas son:

y = (30.5, 32.6, 38.3, 45.7, 33.6, 42.2, 37.5, 38.1, 41.6, 40.4)

Para el estudio de la posible relación lineal entre las variables (x, y), se
desarrollan los siguientes apartados:

Cargar datos y representación de la Nube de Puntos:

plantas=[Link](xi=c(15.3, 17.8, 20.7, 25.1, 16.4, 21.6, 19.6, 18.8, 20.2,

19.4),yi=c(30.5, 32.6, 38.3, 45.7, 33.6, 42.2, 37.5, 38.1, 41.6, 40.4))

plantas
xi yi
1 15.3 30.5
2 17.8 32.6
3 20.7 38.3
4 25.1 45.7
5 16.4 33.6
6 21.6 42.2
7 19.6 37.5
8 18.8 38.1
9 20.2 41.6
10 19.4 40.4

Autor: Francisco Plá Martos 2

 MASTER UNIVERSITARIO IDeA

plot(plantas$xi,plantas$yi)  Figura 1.1

Figura 1.1: Nube de puntos de la muestra (xi, yi) usando el comando plot()

Obtención del coeficiente de correlación y p-valor

Se utiliza el comando [Link](), tal como:

[Link](plantas$xi,plantas$yi,method="pearson",use="[Link]")

Pearson's product-moment correlation

data: plantas$xi and plantas$yi

t = 7.084, df = 8, p-value = 0.0001036
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.7202361 0.9833389
sample estimates:
cor
0.9287109

Autor: Francisco Plá Martos 3

 MASTER UNIVERSITARIO IDeA

Diagramas de dispersión múltiple de todas las variables (en este ejemplo

hay dos variables)

Se procede con el comando pairs([Link]) tal como:

Figura 1.2: Diagrama de dispersión múltiple de la muestra (xi, yi) usando el

comando pairs()

Definición del modelo: comando lm()

Obtención de los coeficientes óptimos y el p-valor asociado a los datos de la

muestra. Consideramos para este ejemplo el Modelo Lineal Simple:

(modeloPlantas=lm(plantas$yi~plantas$xi,data=plantas))

Call:
lm(formula = plantas$yi ~ plantas$xi, data = plantas)

Coefficients:
(Intercept) plantas$xi
7.030 1.592

Por lo tanto, nuestra recta de regresión lineal simple óptima es:

y = 7.03 +1.59 x

Autor: Francisco Plá Martos 4

 MASTER UNIVERSITARIO IDeA

summary(modeloPlantas)

Call:
lm(formula = plantas$yi ~ plantas$xi, data = plantas)

Residuals:
Min 1Q Median 3Q Max
-2.7602 -1.1795 -0.1285 1.0591 2.4932

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 7.0295 4.4182 1.591 0.150262
plantas$xi 1.5916 0.2247 7.084 0.000104 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.857 on 8 degrees of freedom

Multiple R-squared: 0.8625, Adjusted R-squared: 0.8453
F-statistic: 50.18 on 1 and 8 DF, p-value: 0.0001036

Información relevante: p-valor = 0.0001036 < α = 0.05 ( y en realidad para

cualquier α = 0.01, 0.05, 0.10,
Por otro lado F1,N-2 = F1,8 = 50.18 > F1,N-2,α = F1,8,0.05 = 5.0874 (tablas F-Fisher
para α = 0.05) 5.318

Representación de la recta de regresión lineal simple: comando abline()

abline(modeloPlantas,col=c("red"))

Región de confianza al 95% de la media de las imágenes: Es la zona el la

que tenemos una confianza habitualmente al 95% de que la verdadera recta
del modelo de regresión lineal se encuentre.
Es necesario descargar el paquete: visreg

library(visreg)
visreg(lm(yi~xi,data=plantas),type="conditional")

Autor: Francisco Plá Martos 5

 MASTER UNIVERSITARIO IDeA

Información del modelo: comando names(modelo)

names(modeloPlantas)
(residuos_Plantas=modeloPlantas$residuals)
(media_residuos_Plantas=mean(residuos_Plantas))

(predicciones_modeloPlantas=modeloPlantas$[Link])
(coeficientes_modeloPlantas=modeloPlantas$coefficients)
(a=modeloPlantas$coefficients[[1]])
(b=modeloPlantas$coefficients[[2]])

Varianza residual:
(VR=sum((predicciones_modeloPlantas-plantas$yi)^2)/length(plantas$yi))

O equivalentemente:

(VR=sum(residuos_Plantas^2)/ length(plantas$yi))

ANOVA: Comanso anova(modelo)

anova(lm(plantas$yi~plantas$xi,data=plantas))

Analysis of Variance Table

Response: plantas$yi
Df Sum Sq Mean Sq F value Pr(>F)
plantas$xi 1 173.143 173.14 50.184 0.0001036 ***
Residuals 8 27.602 3.45
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Para un nivel de significación α = 0.05, dado que el p-valor = 0.0001036 <

α = 0.05 ( y en realidad para cualquier α = 0.01, 0.05, 0.10, …),
rechazaremos la hipótesis nula H0:{b=0} y por lo tanto existe una relación
lineal entre las variables X e Y.
Por otro lado F1,N-2 = F1,8 = 50.18 > F1,N-2,α = F1,8,0.05 = 5.0874 (tablas F-Fisher
para α = 0.05) , se llega a la misma conclusión: existe una relación lineal
entre las variables X e Y.

Autor: Francisco Plá Martos 6

 MASTER UNIVERSITARIO IDeA

REGRESIÓN LINEAL MÚLTIPLE

Ejemplo2: Concentraciones de ácido acético (Acetic), ácido sulfhídrico

(H2S) y ácido láctico (Lactic) en 30 muestras de queso cheddar maduro.
También se proporciona un valor de sabor subjetivo.
Información sobre los datos:

Uso libre de datos: Mirar más información en la web: [Link]

[Link]

Fuente de datos: David Moore and George McCabe (1989). Introduction to the Practice of Statistics.

El sabor del queso maduro de tipo cheddar está relacionado con la

concentración de diferentes compuestos químicos al final del proceso. Se
considera una tabla de una muestra de 30 datos relacionados con el queso
cheddar en los que se recoge la variable respuesta: gusto (Taste) y tres
variables explicativas relacionadas con la concentración de ácido acético
(Acetic), ácido sulfhídrico (H2S) y ácido láctico (Lactic)

Cargar datos y representación de la Nube de Puntos:

Queso=[Link]("[Link]

head(queso)

Case Acetic H2S Lactic Taste

1 1 4.543 3.135 0.86 12.3
2 2 5.159 5.043 1.53 20.9
3 3 5.366 5.438 1.57 39.0
4 4 5.759 7.496 1.81 47.9
5 5 4.663 3.807 0.99 5.6
6 6 5.697 7.601 1.09 25.9

Autor: Francisco Plá Martos 7

 MASTER UNIVERSITARIO IDeA

pairs(queso)  Figura 2.1

Figura 2.1: Diagrama de dispersión múltiple usando el comando pairs()

Autor: Francisco Plá Martos 8

 MASTER UNIVERSITARIO IDeA

Obtención del coeficiente de correlación de todas las combinaciones de

variables:
cor(queso)
Case Acetic H2S Lactic Taste
Case 1.00000000 0.2838356 0.04383232 0.05658513 -0.2148631
Acetic 0.28383559 1.0000000 0.61795591 0.60434768 0.5495295
H2S 0.04383232 0.6179559 1.00000000 0.64389703 0.7557630
Lactic 0.05658513 0.6043477 0.64389703 1.00000000 0.7034822
Taste -0.21486308 0.5495295 0.75576301 0.70348216 1.0000000

Instalación del paquete tree

[Link]("tree")
library(tree)
modeloT=tree(queso$Taste~., data=queso)
plot(modeloT)
text(modeloT)

Figura 2.2: Diagrama de árbol del modelo

Autor: Francisco Plá Martos 9

 MASTER UNIVERSITARIO IDeA

Con la información proporcionada por el árbol de decisión podemos

construir el siguiente modelo:

Definición del modelo: comando lm()

(Modelo1=lm(queso$Taste~queso$H2S+queso$Lactic, data=queso))

Call:
lm(formula = queso$Taste ~ queso$H2S + queso$Lactic, data = queso)

Coefficients:
(Intercept) queso$H2S queso$Lactic
-27.622 3.953 19.884

Por lo tanto, nuestra recta de regresión lineal múltiple óptima es:

y = -27.622 +3.953 x1 + 19.884 x2

es decir:

Taste= -27.622 +3.953 H2 S + 19.884 Lactic

summary(Modelo1)

Call:
lm(formula = queso$Taste ~ queso$H2S + queso$Lactic, data = queso)

Residuals:
Min 1Q Median 3Q Max
-17.363 -6.547 -1.161 4.848 25.606

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -27.622 9.001 -3.069 0.00485 **
queso$H2S 3.953 1.135 3.483 0.00170 **
queso$Lactic 19.884 7.970 2.495 0.01902 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 9.946 on 27 degrees of freedom

Multiple R-squared: 0.6515, Adjusted R-squared: 0.6257
F-statistic: 25.24 on 2 and 27 DF, p-value: 6.602e-07

Autor: Francisco Plá Martos 10

 MASTER UNIVERSITARIO IDeA

Para un nivel de significación α = 0.05, el p-valor = 6.602e-07 es muy

significativo, por lo que rechazaremos la hipótesis nula H0:{b1 = b2 = 0} y
por lo tanto existe una relación lineal entre las variable respuesta Y (Taste)
y las variables explicativas, H2S y ácido láctico (Lactic). Además se puede
observar que la variable H2S es más significativa (p-valor=0.00170) que el
ácido láctico (p-valor = 0.01902). También obtenemos la F-Fisher de la tabla
ANOVA asociada a la regresión lineal múltiple con k=2: Fk,N-(k+1) = F2,27 =
25.24 . El valor R2 = 0.6515

Es muy importante estudiar el error de los residuos 𝜀𝜀𝑖𝑖 = 𝑦𝑦𝑖𝑖 − 𝑦𝑦�𝚤𝚤 , que tiene
que seguir una distribución normal de media 0 y cumple homocedasticidad:
𝜎𝜎12 = 𝜎𝜎22 = ⋯ = 𝜎𝜎𝑁𝑁2

Para ello, en primer lugar estudiamos si nuestros residuos siguen un

distribución Normal:

[Link](Modelo1$residuals)

Shapiro-Wilk normality test

data: Modelo1$residuals
W = 0.97947, p-value = 0.8114

Dado que p-valor = 0.8114 > α = 0.05, entonces no podemos rechazar la

hipótesis de normalidad

Con el comando estudiamos el resto de las condiciones del error:

plot(Modelo1)
Hit <Return> to see next plot: y
Hit <Return> to see next plot: y
Hit <Return> to see next plot: y
Hit <Return> to see next plot: y

Autor: Francisco Plá Martos 11

 MASTER UNIVERSITARIO IDeA

(a) (b)

(c) (d)

(a) Los residuos se deben de comportar de forma aleatoria, sin ningún patrón especifico
(b) Normalidad de los residuos, los datos están alineados sobre la recta, muy pocos se
escapan
(c) Homocedasticidad. Se debe vigilar la estructura.
(d) Detecta valores atípicos: aquellos que están fuera de la región punteada, los que tienen
distancia de Cook mayor a 1, D>1

Autor: Francisco Plá Martos 12

 MASTER UNIVERSITARIO IDeA

Construcción de la tabla ANOVA:

La Tabla ANOVA del modelo de regresión lineal múltiple nos indica la
formulación que tenemos que calcular

DF_Explicada=2
(DF_Residual=length(queso$Case)-(DF_Explicada+1))
[1] 27
(Explicada=sum((Modelo1$[Link]-mean(queso$Taste))^2))
[1] 4992.828
(Residual=sum((Modelo1$[Link]-queso$Taste)^2))
#(Residual=sum((Modelo1$residuals)^2))
[1] 2670.702
(A=Explicada/DF_Explicada)
[1] 2496.414
(B=Residual/DF_Residual)
[1] 98.91489
(F=A/B)
[1] 25.238

Por lo tanto:

Fuente de Suma de GL Media de F p-value

variación cuadrados cuadrados
Explicativa 4992.828 2 2496.414 25.238 6.602e-07
Regresora
No 2670.702 30 - (2+1) 98.91489
explicativa = 27
Residual

Autor: Francisco Plá Martos 13

 MASTER UNIVERSITARIO IDeA

T7.3 Estimación y contraste de

hipótesis
ANOVA de 1 factor

Ejemplo:

Dada un archivo de datos teórico basado en 100 consumidores de

yogures se les ha pedido su opinión de percepción sobre los siguientes
descriptores: Color, Sabor, Dulce, Acido, Olor y Cremoso en 5 muestras
de marcas diferentes de yogur de plátano (mirar archivo de diseño de
registro). Se pretende analizar si existen diferencias significativas
respecto de los descriptores entre las 5 marcas de yogures

Observación: Tenemos un archivo de datos .txt que se tiene que

descargar en R y que tiene un tamaño de 5 × 100 = 500 filas.

Cargamos los datos y creamos la tabla de datos:

setwd("…")

Yogur=[Link]("Datos_yogures_5_muestras_100_consum_mIDeA201920.txt",header = T)

MUESTRA=factor(Yogur$MUESTRA, labels=c("Muestra 1","Muestra 2","Muestra

3","Muestra 4","Muestra 5"))

Vamos a estudiar la variable Color:

Representamos el diagrama QQ para verificar si los datos recogidos de la

variable Color siguen una distribución Normal:

qqnorm(Yogur$Color)
qqline(Yogur$Color,col="red")

Autor: Francisco Plá Martos 14

 MASTER UNIVERSITARIO IDeA

Observamos que hay muchos datos de la variable Color que se alejan de la

línea diagonal, por lo que no parece que sigan una distribución Normal

Test de Normalidad:

# Normalidad
[Link](Yogur$Color) # p-value = 2.2e-16 -> rechazamos normalidad

Shapiro-Wilk normality test

data: Yogur$Color
W = 0.87655, p-value < 2.2e-16

Rechazamos normalidad dado que α = 0.05 > p-valor = 2.2e-16

# homecedasticidad:

library(carData)
library(car)
leveneTest(y=Yogur$Color,group=MUESTRA,center="median
Levene's Test for Homogeneity of Variance (center = "median")
Df F value Pr(>F)
group 4 4.3466 0.00183 **
495

Autor: Francisco Plá Martos 15

 MASTER UNIVERSITARIO IDeA

Rechazamos homocedasticidad dado que α = 0.05 > p-valor = 0.00183

Por lo tanto consideramos un test no paramétrico alternativo a la ANOVA,

describimos a continuación el Test de tipo Kruskal-Wallis:

[Link](Yogur$Color~MUESTRA,data=Yogur)

Kruskal-Wallis rank sum test

data: Yogur$Color by MUESTRA

Kruskal-Wallis chi-squared = 154.96, df = 4, p-value < 2.2e-16

Dado que el p-valor es muy significativo es decir, muy pequeño, entonces se

rechaza la hipótesis de igualdad medias entre las 5 muestras de marcas
diferentes, por lo tanto sí que existen diferencias significativas entre las
marcas respecto de la variable Color

Con el siguiente comando podemos observar la variabilidad de los datos de

Color para cada una de las 5 marcas:
library(ggplot2)
ggplot(data=Yogur, mapping = aes(x = Yogur$Color, colour = MUESTRA))
+ geom_histogram()+theme_bw() + facet_grid(. ~ MUESTRA) +
theme([Link] = "none")

Autor: Francisco Plá Martos 16

 MASTER UNIVERSITARIO IDeA

Con un diagrama de cajas comparamos las 5 muestras de marcas respecto

de la variable Color

boxplot(Yogur$Color~MUESTRA,[Link]=T,col=c("red","cyan","pink"
,"brown","orange"))

Se observan que en efecto existen diferencias significativas entre las

muestras

POST-HOC
En el caso de que nuestros datos siguieran una distribución Normal,
realizaríamos el test Post-Hoc
[Link](Yogur$Color,MUESTRA,[Link]="bonferroni")

Autor: Francisco Plá Martos 17

 MASTER UNIVERSITARIO IDeA

Pero dado que no es así, realizaremos un test Post-Hoc no paramétrico tal

como el siguiente:

[Link](Yogur$Color,MUESTRA,[Link] = "holm")

Pairwise comparisons using Wilcoxon rank sum test

data: Yogur$Color and MUESTRA

Muestra 1 Muestra 2 Muestra 3 Muestra 4

Muestra 2 1.1e-11 - - -
Muestra 3 0.21 < 2e-16 - -
Muestra 4 6.0e-11 0.91 < 2e-16 -
Muestra 5 3.0e-08 0.28 7.0e-16 0.28

P value adjustment method: holm

El post-hoc nos indica que existen diferencias significativas entre las

muestras 1-2, 1-4, 1-5, 2-3, 3-4 y 3-5 a un nivel de confianza del 95%

Para estudiar la normalidad de todas las variables descriptivas al mismo

tiempo, R te ofrece la posibilidad de hacer un bucle “for” tal y como se indica
a continuación:
for(i in 3:8)
{print([Link](Yogur[,i]))}

COLOR
Shapiro-Wilk normality test

data: Yogur[, i]
W = 0.87655, p-value < 2.2e-16

SABOR
Shapiro-Wilk normality test

data: Yogur[, i]
W = 0.87979, p-value < 2.2e-16

DULCE
Shapiro-Wilk normality test

data: Yogur[, i]
W = 0.90712, p-value < 2.2e-16

Autor: Francisco Plá Martos 18

 MASTER UNIVERSITARIO IDeA

ACIDO
Shapiro-Wilk normality test

data: Yogur[, i]
W = 0.91058, p-value < 2.2e-16

OLOR
Shapiro-Wilk normality test

data: Yogur[, i]
W = 0.90005, p-value < 2.2e-16

CREMOSO
Shapiro-Wilk normality test

data: Yogur[, i]
W = 0.9149, p-value = 3.709e-16

En todos los casos no siguen una distribución Normal, por lo que se tendrá
que realizar un test no paramétrico de tipo Kruskal-Wallis. Proponemos
realizar otro bucle:
for(i in 3:8)
{print([Link](Yogur[,i]~MUESTRA,data=Yogur))}

COLOR
Kruskal-Wallis rank sum test

data: Yogur[, i] by MUESTRA

Kruskal-Wallis chi-squared = 154.96, df = 4, p-value < 2.2e-16

SABOR
Kruskal-Wallis rank sum test

data: Yogur[, i] by MUESTRA

Kruskal-Wallis chi-squared = 76.877, df = 4, p-value = 7.985e-16

DULCE
Kruskal-Wallis rank sum test

data: Yogur[, i] by MUESTRA

Kruskal-Wallis chi-squared = 48.298, df = 4, p-value = 8.178e-10

Autor: Francisco Plá Martos 19

 MASTER UNIVERSITARIO IDeA

ACIDO
Kruskal-Wallis rank sum test

data: Yogur[, i] by MUESTRA

Kruskal-Wallis chi-squared = 146.4, df = 4, p-value < 2.2e-16

OLOR
Kruskal-Wallis rank sum test

data: Yogur[, i] by MUESTRA

Kruskal-Wallis chi-squared = 134.5, df = 4, p-value < 2.2e-16

CREMOSO
Kruskal-Wallis rank sum test

data: Yogur[, i] by MUESTRA

Kruskal-Wallis chi-squared = 3.8497, df = 4, p-value = 0.4267

Todos los descriptivos son significativo a excepción de Cremoso. Si

representamos Cremoso en diagramas de cajas vemos claro el resultado de
dicho test (una imagen vale más que mil palabras)

Autor: Francisco Plá Martos 20

 MASTER UNIVERSITARIO IDeA

T7.4 ANOVA de 2 factores

ANOVA de 2 factores con interacciones

Ejemplo1:

Medimos el rendimiento en Kg/ha de una producción agrícola teniendo

en cuenta el tipo de semilla y el fertilizante. Consideramos dos tipos de
semilla y tres tipos de fertilizante y se recoge el rendimiento en la tabla
de a continuación en el que se toma una muestra de tamaño 2 por semilla
y fertilizante.

Semilla 1 Semilla 2
Fertilizante 1 (14.3, 14.5 ) (12.6, 11.2 )
Fertilizante 2 (18.1, 17.6 ) (10.5, 12.8 )
Fertilizante 3 (17.6, 18.2 ) (15.7, 17.5 )

Tenemos por lo tanto un diseño ANOVA de 2 factores con interacciones, tal

que:
• Factor Fertilizante con 3 niveles
• Factor Semilla con 2 niveles
• Tamaño de la muestra en cada cruce k = 2
• Tenemos 12 datos

Procedemos a cargar los datos con R:

rendimiento=c(14.3,14.5,18.1,17.6,17.6,18.2,12.6,11.2,10.5,12.8,15.7,17.5)

semilla=c(1,1,1,1,1,1,2,2,2,2,2,2)

fertilizante=c(1,1,2,2,3,3,1,1,2,2,3,3)

fertilizante=
factor(fertilizante,labels=c("Fertilizante 1","Fertilizante 2","Fertilizante
3"))

semilla=factor(semilla,labels=c("Semilla 1","Semilla 2"))

Autor: Francisco Plá Martos 21

 MASTER UNIVERSITARIO IDeA

(datos=[Link](rendimiento,fertilizante,semilla))

rendimiento fertilizante semilla

1 14.3 Fertilizante 1 Semilla 1
2 14.5 Fertilizante 1 Semilla 1
3 18.1 Fertilizante 2 Semilla 1
4 17.6 Fertilizante 2 Semilla 1
5 17.6 Fertilizante 3 Semilla 1
6 18.2 Fertilizante 3 Semilla 1
7 12.6 Fertilizante 1 Semilla 2
8 11.2 Fertilizante 1 Semilla 2
9 10.5 Fertilizante 2 Semilla 2
10 12.8 Fertilizante 2 Semilla 2
11 15.7 Fertilizante 3 Semilla 2
12 17.5 Fertilizante 3 Semilla 2

Con el siguiente comando obtenemos una tabla de medias en cada una de las
celdas de factores

tapply(datos$rendimiento,list(datos$fertilizante,datos$semilla),mean)
Semilla 1 Semilla 2
Fertilizante 1 14.40 11.90
Fertilizante 2 17.85 11.65
Fertilizante 3 17.90 16.60

Antes de aplicar una tabla ANOVA de 2 factores con interacciones tenemos

que estudiar si la variable cuantitativa rendimiento sigue una distribución
normal:

qqnorm(datos$rendimiento)
qqline(datos$rendimiento)

En la Figura 4.1, se observa que los datos de la variable Rendimiento están

alineados sobre la recta, muy pocos se escapan, por lo que parece que en
efecto siguen una distribución normal aunque seguramente el p-valor que
obtengamos no sea muy grande respecto del α = 0.05

Autor: Francisco Plá Martos 22

 MASTER UNIVERSITARIO IDeA

Con el siguiente test de normalidad, obtenemos el correspondiente p-valor:

[Link](datos$rendimiento)

Shapiro-Wilk normality test

data: datos$rendimiento
W = 0.89471, p-value = 0.1355

Figura 4.1: Diagrama QQ de Normalidad de la variable Rendimiento.

ANOVA dos factores con interacciones

En efecto, dado que α = 0.05 < p-value = 0.1355, no podemos rechazar que
los datos de la variable rendimiento sigan una distribución Normal.
Procedemos a obtener la tabla ANOVA de 2 factores con interacciones de
nuestro experimento:

modelo=lm(datos$rendimiento~semilla*fertilizante)

Autor: Francisco Plá Martos 23

 MASTER UNIVERSITARIO IDeA

# o equivalentemente:
#modelo=lm(datos$rendimiento~semilla+fertilizante+semilla*fertilizante)
anova(modelo)

Analysis of Variance Table

Response: datos$rendimiento
Df Sum Sq Mean Sq F value Pr(>F)
semilla 1 33.333 33.333 35.9066 0.0009711 ***
fertilizante 2 34.160 17.080 18.3986 0.0027556 **
semilla:fertilizante 2 13.047 6.523 7.0269 0.0267830 *
Residuals 6 5.570 0.928
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Por lo tanto, todos los factores son significativos para un 95 % de confianza,

aunque el que más es el factor Semilla y las interacciones las que menos
significativas.

ANOVA de 2 factores sin interacciones

Ejemplo2:

Se considera el mismo experimento que en Ejemplo 1 excepto que

tomamos sólo una medición por celda. Tenemos la siguiente tabla:

Semilla 1 Semilla 2
Fertilizante 1 14.3 12.6
Fertilizante 2 18.1 10.5
Fertilizante 3 17.6 15.7

Tenemos por lo tanto un diseño ANOVA de 2 factores, tal que:

• Factor Fertilizante con 3 niveles
• Factor Semilla con 2 niveles

• Tenemos 6 datos

Procedemos a cargar los datos con R:

rendimiento=c(14.3,18.1,17.6,12.6,10.5,15.7)

semilla=c(1,1,1,2,2,2)

Autor: Francisco Plá Martos 24

 MASTER UNIVERSITARIO IDeA

fertilizante=c(1,2,3,1,2,3)
fertilizante=
factor(fertilizante,labels=c("Fertilizante 1","Fertilizante 2","Fertilizante
3"))

semilla=factor(semilla,labels=c("Semilla 1","Semilla 2"))

(datos=[Link](rendimiento,fertilizante,semilla))

rendimiento fertilizante semilla

1 14.3 Fertilizante 1 Semilla 1
2 18.1 Fertilizante 2 Semilla 1
3 17.6 Fertilizante 3 Semilla 1
4 12.6 Fertilizante 1 Semilla 2
5 10.5 Fertilizante 2 Semilla 2
6 15.7 Fertilizante 3 Semilla 2

tapply(datos$rendimiento,list(datos$fertilizante,datos$semilla),mean)
Semilla 1 Semilla 2
Fertilizante 1 14.3 12.6
Fertilizante 2 18.1 10.5
Fertilizante 3 17.6 15.7

Estudiamos Normalidad:

Diagrama

qqnorm(datos$rendimiento)
qqline(datos$rendimiento)

En la Figura 4.2, se que en efecto la variable Rendiemiento sigue una

distribución normal.

Test de normalidad:
[Link](datos$rendimiento)

Shapiro-Wilk normality test

data: datos$rendimiento
W = 0.95398, p-value = 0.7723

Autor: Francisco Plá Martos 25

 MASTER UNIVERSITARIO IDeA

Figura 4.2: Diagrama QQ de Normalidad de la variable Rendimiento.

ANOVA dos factores

Por lo tanto, no podemos rechazar normalidad de la variable Rendimiento,

dado que α = 0.05 < p-value = 0.7723.

Procedemos a obtener la tabla ANOVA de 2 factores:

modelo=lm(datos$rendimiento~semilla+fertilizante)
anova(modelo)

Response: datos$rendimiento
Df Sum Sq Mean Sq F value Pr(>F)
semilla 1 20.907 20.9067 3.7256 0.1933
fertilizante 2 10.990 5.4950 0.9792 0.5053
Residuals 2 11.223 5.6117

Es muy interesante observar por medio de los resultados de la Tabla

ANOVA de 2 factores, que ninguno de los factores son significativos en este
experimento particular, si no se consideran interacciones, es decir, si no se
replica el experimento.

Autor: Francisco Plá Martos 26