¿Han visto una campana? De esas que hace 30 años eran repicadas con una piedra o con un mazo, para disciplinarnos a ingresar a clases. Una distribución normal en estadística debe parecerse a esas campanas, vista en dos dimensiones.
 |
|
Si al graficar los datos obtenidos de un experimento, la figura carece de esa silueta, hay que transformar los valores.
 |
valores
sin transformar (%) <-------> volares transformados (grados)
Transformación
angular
Excel: =DEGREES(ASIN(SQRT(D2/100)))
SAS : ARCSEN = ARSIN(SQRT(PHA/100))*180/constant("pi") R : ARCSEN <- asin(sqrt(PHA/100))*(180/pi) |
Algo parecido a la trasformación de Optimus Prime, de robot a tráiler, en las películas de Transformers. El robot, datos no normales, cambia a tráiler, datos normales, para poder rodar o, en este caso, correr análisis paramétricos de ANOVA.
 |
|
Encontraran
los fundamentos del análisis de varianza del Diseño de Bloques Divididos, y como ajustar datos porcentuales
a una distribución normal en el Capítulo 12 de Experimentación en Agricultura. La siguiente es nuestra propuesta de códigos.
Códigos
en R
#***********************************************************
# EXPERIMENTACIÓN EN AGRICULTURA
# Fernández Escobar, R.; Trapero, A.; Domínguez, J. 2010
# CAPÍTULO 12: DISEÑO EN BLOQUES DIVIDIDOS
#************************************************************
# Importar datos
hojas <- read.csv("C:/Users/Administrator/Desktop/tabla 12.1.csv", header=TRUE)
attach(hojas); hojas## SUELO TRAT REP PHA
## 1 1 1 1 85.5
## 2 1 1 2 87.8
## 3 1 1 3 81.4
## 4 1 1 4 74.3
## 5 1 1 5 94.5
## 6 2 1 1 89.3
## 7 2 1 2 76.8
## 8 2 1 3 70.0
## 9 2 1 4 78.5
## 10 2 1 5 81.4
## 11 1 2 1 16.7
## 12 1 2 2 12.0
## 13 1 2 3 38.1
## 14 1 2 4 22.2
## 15 1 2 5 18.0
## 16 2 2 1 6.9
## 17 2 2 2 24.8
## 18 2 2 3 10.0
## 19 2 2 4 14.1
## 20 2 2 5 9.2
## 21 1 3 1 65.5
## 22 1 3 2 60.7
## 23 1 3 3 38.2
## 24 1 3 4 78.7
## 25 1 3 5 56.9
## 26 2 3 1 52.3
## 27 2 3 2 68.7
## 28 2 3 3 41.5
## 29 2 3 4 28.1
## 30 2 3 5 38.6
## 31 1 4 1 6.2
## 32 1 4 2 7.7
## 33 1 4 3 5.4
## 34 1 4 4 3.6
## 35 1 4 5 15.0
## 36 2 4 1 0.7
## 37 2 4 2 6.3
## 38 2 4 3 11.2
## 39 2 4 4 2.9
## 40 2 4 5 4.5#### Etiquetar y transformar datos
# Etiquetar variables categóricas
SUELO <- factor(SUELO, levels = c(1,2), labels = c("Cubierta Vegetal", "Suelo Desnudo"))
TRAT <- factor(TRAT, levels = c(1,2,3,4), labels = c("Sin tratar", "Primavera", "Otoño", "Prim + Otoño" ))
hojasID <- data.frame(SUELO, TRAT, REP, PHA); hojasID# Porcentajes transformados a grados
help(arcsin)ARCSEN <- asin(sqrt(PHA/100))*(180/pi)
hojasTRANS <- data.frame(hojasID, ARCSEN); hojasTRANS## SUELO TRAT REP PHA ARCSEN
## 1 Cubierta Vegetal Sin tratar 1 85.5 67.617458
## 2 Cubierta Vegetal Sin tratar 2 87.8 69.556412
## 3 Cubierta Vegetal Sin tratar 3 81.4 64.451360
## 4 Cubierta Vegetal Sin tratar 4 74.3 59.539005
## 5 Cubierta Vegetal Sin tratar 5 94.5 76.436623
## 6 Suelo Desnudo Sin tratar 1 89.3 70.906630
## 7 Suelo Desnudo Sin tratar 2 76.8 61.205877
## 8 Suelo Desnudo Sin tratar 3 70.0 56.789089
## 9 Suelo Desnudo Sin tratar 4 78.5 62.375113
## 10 Suelo Desnudo Sin tratar 5 81.4 64.451360
## 11 Cubierta Vegetal Primavera 1 16.7 24.120456
## 12 Cubierta Vegetal Primavera 2 12.0 20.267901
## 13 Cubierta Vegetal Primavera 3 38.1 38.115736
## 14 Cubierta Vegetal Primavera 4 22.2 28.110190
## 15 Cubierta Vegetal Primavera 5 18.0 25.104090
## 16 Suelo Desnudo Primavera 1 6.9 15.229055
## 17 Suelo Desnudo Primavera 2 24.8 29.867504
## 18 Suelo Desnudo Primavera 3 10.0 18.434949
## 19 Suelo Desnudo Primavera 4 14.1 22.055199
## 20 Suelo Desnudo Primavera 5 9.2 17.656820
## 21 Cubierta Vegetal Otoño 1 65.5 54.029615
## 22 Cubierta Vegetal Otoño 2 60.7 51.178433
## 23 Cubierta Vegetal Otoño 3 38.2 38.174712
## 24 Cubierta Vegetal Otoño 4 78.7 62.514816
## 25 Cubierta Vegetal Otoño 5 56.9 48.966066
## 26 Suelo Desnudo Otoño 1 52.3 46.318268
## 27 Suelo Desnudo Otoño 2 68.7 55.981260
## 28 Suelo Desnudo Otoño 3 41.5 40.106090
## 29 Suelo Desnudo Otoño 4 28.1 32.011829
## 30 Suelo Desnudo Otoño 5 38.6 38.410324
## 31 Cubierta Vegetal Prim + Otoño 1 6.2 14.418226
## 32 Cubierta Vegetal Prim + Otoño 2 7.7 16.110382
## 33 Cubierta Vegetal Prim + Otoño 3 5.4 13.437174
## 34 Cubierta Vegetal Prim + Otoño 4 3.6 10.937416
## 35 Cubierta Vegetal Prim + Otoño 5 15.0 22.786498
## 36 Suelo Desnudo Prim + Otoño 1 0.7 4.799319
## 37 Suelo Desnudo Prim + Otoño 2 6.3 14.536577
## 38 Suelo Desnudo Prim + Otoño 3 11.2 19.552108
## 39 Suelo Desnudo Prim + Otoño 4 2.9 9.804905
## 40 Suelo Desnudo Prim + Otoño 5 4.5 12.247324attach(hojasTRANS)#### Histogramas y gráficos de densidades
# Datos sin transformar
par(mfrow=c(1,2))
hist(PHA, prob=TRUE, col="grey")# prob=TRUE para probabilidades, no recuentos
lines(density(PHA), col="blue", lwd=2)
lines(density(PHA, adjust=2), lty="dotted", col="darkgreen", lwd=2)
# Datos transformados
hist(ARCSEN, prob=TRUE, col="grey")
lines(density(ARCSEN), col="blue", lwd=2)
lines(density(ARCSEN, adjust=2), lty="dotted", col="darkgreen", lwd=2)
#### Gráficos Box and Whisker
library(HH)bwplot(PHA ~ TRAT, data=hojasTRANS,
scales=list(cex=1.5),
ylab=list("Hojas enfermas (%)", cex=1.5),
xlab=list("Momentos de los tratamientos",cex=1.5))bwplot(ARCSEN ~ TRAT, data=hojasTRANS,
scales=list(cex=1.5),
ylab=list("ARCSEN, datos transformados", cex=1.5),
xlab=list("Momentos de los tratamientos",cex=1.5))
#### Estadísticas descriptivas por tratamientos
library(plyr)source("C:/Users/Administrator/Desktop/summarySE.r") # Función summarySE
# Hojas enfermas (%)
summarySE(hojasTRANS, measurevar="PHA", groupvars=c("SUELO", "TRAT"))## SUELO TRAT N PHA sd se ci
## 1 Cubierta Vegetal Sin tratar 5 84.70 7.505665 3.356635 9.319513
## 2 Cubierta Vegetal Primavera 5 21.40 10.019232 4.480737 12.440519
## 3 Cubierta Vegetal Otoño 5 60.00 14.707821 6.577538 18.262172
## 4 Cubierta Vegetal Prim + Otoño 5 7.58 4.403635 1.969365 5.467835
## 5 Suelo Desnudo Sin tratar 5 79.20 7.031003 3.144360 8.730143
## 6 Suelo Desnudo Primavera 5 13.00 7.090487 3.170962 8.804002
## 7 Suelo Desnudo Otoño 5 45.84 15.413241 6.893011 19.138067
## 8 Suelo Desnudo Prim + Otoño 5 5.12 3.975173 1.777751 4.935829summarySE(hojasTRANS, measurevar="PHA", groupvars=c("SUELO", "SUELO"), na.rm=TRUE)## SUELO N PHA sd se ci
## 1 Cubierta Vegetal 20 43.42 32.69583 7.31101 15.30212
## 2 Suelo Desnudo 20 35.79 31.32146 7.00369 14.65889# ARCSEN, datos transformados
summarySE(hojasTRANS, measurevar="ARCSEN", groupvars=c("SUELO", "TRAT"))## SUELO TRAT N ARCSEN sd se ci
## 1 Cubierta Vegetal Sin tratar 5 67.52017 6.260501 2.799781 7.773438
## 2 Cubierta Vegetal Primavera 5 27.14367 6.743447 3.015761 8.373095
## 3 Cubierta Vegetal Otoño 5 50.97273 8.809407 3.939686 10.938323
## 4 Cubierta Vegetal Prim + Otoño 5 15.53794 4.463773 1.996260 5.542506
## 5 Suelo Desnudo Sin tratar 5 63.14561 5.164785 2.309762 6.412927
## 6 Suelo Desnudo Primavera 5 20.64871 5.705023 2.551364 7.083722
## 7 Suelo Desnudo Otoño 5 42.56555 9.065871 4.054381 11.256765
## 8 Suelo Desnudo Prim + Otoño 5 12.18805 5.477709 2.449706 6.801473summarySE(hojasTRANS, measurevar="ARCSEN", groupvars=c("SUELO", "SUELO"), na.rm=TRUE)## SUELO N ARCSEN sd se ci
## 1 Cubierta Vegetal 20 40.29363 21.68466 4.848836 10.148731
## 2 Suelo Desnudo 20 34.63698 21.22787 4.746696 9.934948#### Gráficos de barras para Hojas enfermas (%)
library(ggplot2)gPHA <- summarySE(hojasTRANS, measurevar="PHA", groupvars=c("SUELO", "TRAT"))
# Barras de error representan errores típicos de la media
ggplot(gPHA, aes(x = TRAT, y = PHA, fill = SUELO)) +
geom_bar(position = position_dodge(), stat = "identity") +
geom_errorbar(aes(ymin = PHA-se, ymax = PHA+se),
width = .2, # Ancho de las barras de error
position = position_dodge(.9))
# Gráfico de %, pero más elaborado
ggplot(gPHA, aes(x = TRAT, y = PHA, fill = SUELO)) +
geom_bar(position = position_dodge(), stat ="identity",
colour = "black", # utilice esquemas negros,
size = .3) + # líneas más finas
geom_errorbar(aes(ymin = PHA-se, ymax=PHA+se),
size = .3, # líneas más finas
width = .2,
position = position_dodge(.9)) +
xlab("Momentos de los tratamientos") +
ylab("Hojas enfermas (%)") +
scale_fill_hue(name = "SUELO", # etiqueta de leyenda, use colores oscuros
breaks = c("Cubierta Vegetal", "Suelo Desnudo"),
labels = c("Vegetal", "Desnudo")) +
ggtitle("Porcentaje de hojas enfermas \nAplicación de tratamientos y dos sistemas de manejo") +
scale_y_continuous(breaks=0:20*4) +
theme_bw()
#### Gráfico de barras para Hojas enfermas (ARCSEN)
gARCSEN <- summarySE(hojasTRANS, measurevar="ARCSEN", groupvars=c("SUELO", "TRAT"))
# Barras de error representan errores típicos de la media
ggplot(gARCSEN, aes(x = TRAT, y = ARCSEN, fill = SUELO)) +
geom_bar(position = position_dodge(), stat = "identity") +
geom_errorbar(aes(ymin = ARCSEN-se, ymax = ARCSEN+se),
width=.2, # Ancho de las barras de error
position=position_dodge(.9))
# Gráfico de ARCSEN, pero más elaborado
ggplot(gARCSEN, aes(x = TRAT, y = ARCSEN, fill = SUELO)) +
geom_bar(position = position_dodge(), stat = "identity",
colour = "black", # utilice esquemas negros,
size = .3) + # líneas más finas
geom_errorbar(aes(ymin = ARCSEN-se, ymax=ARCSEN+se),
size = .3, # líneas más finas
width = .2,
position=position_dodge(.9)) +
xlab("Momentos de los tratamientos") +
ylab("Hojas enfermas (ARCSEN)") +
scale_fill_hue(name = "SUELO", # etiqueta de leyenda, use colores oscuros
breaks = c("Cubierta Vegetal", "Suelo Desnudo"),
labels = c("Vegetal", "Desnudo")) +
ggtitle("ARCSEN de hojas enfermas \nAplicación de tratamientos y dos sistemas de manejo") +
scale_y_continuous(breaks=0:20*4) + theme_bw()
#### Diversos gráfico para comparar tratamientos
interaction2wt(ARCSEN ~ SUELO + TRAT + REP, data = hojasTRANS,
par.strip.text=list(cex=.7))
interaction2wt(ARCSEN ~ SUELO + TRAT, data = hojasTRANS)
#### Prueba de Bartlett
bartlett.test(PHA ~ TRAT, data = hojasTRANS) # para PHA## Bartlett's K-squared = 14.675, df = 3, p-value = 0.002116bartlett.test(ARCSEN ~ TRAT, data = hojasTRANS) # para ARCSEN## Bartlett's K-squared = 4.0223, df = 3, p-value = 0.2591
#### Análisis de Varianza para Diseño en Bloques Divididos
#### ANOVA para PHA y ARCSEN (opción 1 )
library(agricolae)modelo1.aov<-with(data = hojasTRANS, strip.plot(REP, SUELO, TRAT, PHA)) # sin transfor## Response: PHA
## Df Sum Sq Mean Sq F value Pr(>F)
## REP 4 184 46.1 0.3962 0.80765
## SUELO 1 582 582.2 7.2630 0.05438 .
## Ea 4 321 80.2 0.6895 0.61309
## TRAT 3 35783 11927.6 132.2627 1.853e-09 ***
## Eb 12 1082 90.2 0.7757 0.66649
## TRAT:SUELO 3 186 62.1 0.5340 0.66764
## Ec 12 1395 116.3
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## cv(a) = 22.6 %, cv(b) = 24 %, cv(c) = 27.2 %, Mean = 39.605modelo2.aov<-with(data = hojasTRANS, strip.plot(REP, SUELO, TRAT, ARCSEN)) # transfor## Response: ARCSEN
## Df Sum Sq Mean Sq F value Pr(>F)
## REP 4 84.1 21.0 0.3896 0.81212
## SUELO 1 320.0 320.0 8.4218 0.04403 *
## Ea 4 152.0 38.0 0.7045 0.60394
## TRAT 3 16043.5 5347.8 120.8040 3.135e-09 ***
## Eb 12 531.2 44.3 0.8208 0.63109
## TRAT:SUELO 3 38.1 12.7 0.2354 0.87000
## Ec 12 647.2 53.9
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## cv(a) = 16.5 %, cv(b) = 17.8 %, cv(c) = 19.6 %, Mean = 37.4653#### Comparación de medias LSD
# Cobertura de suelo
comparación1<-LSD.test(PHA, SUELO, modelo2.aov$gl.a, modelo2.aov$Ea); comparación1## $groups
## trt means M
## 1 Cubierta Vegetal 43.42 a
## 2 Suelo Desnudo 35.79 bcomparación1<-LSD.test(ARCSEN, SUELO, modelo2.aov$gl.a, modelo2.aov$Ea); comparación1## $groups
## trt means M
## 1 Cubierta Vegetal 40.29363 a
## 2 Suelo Desnudo 34.63698 b# Momentos de los tratamientos
comparación2<-LSD.test(PHA, TRAT, modelo2.aov$gl.b, modelo2.aov$Eb); comparación2## $groups
## trt means M
## 1 Sin tratar 81.95 a
## 2 Otoño 52.92 b
## 3 Primavera 17.20 c
## 4 Prim + Otoño 6.35 dcomparación2<-LSD.test(ARCSEN, TRAT, modelo2.aov$gl.b, modelo2.aov$Eb); comparación2## $groups
## trt means M
## 1 Sin tratar 65.33289 a
## 2 Otoño 46.76914 b
## 3 Primavera 23.89619 c
## 4 Prim + Otoño 13.86299 d#### ANOVA para PHA (opción 2)
modelo3.aov <- aov(PHA ~ SUELO * TRAT + REP:SUELO + REP:TRAT, data = hojasTRANS) # sin transformar
summary(modelo3.aov)shapiro.test(modelo3.aov$residuals)# Gráfico de diagnósticos
par(mfrow=c(1,2))
plot(modelo3.aov, 1)
plot(modelo3.aov, 2)
plot(modelo3.aov, 3)
plot(modelo3.aov, 4)
plot(modelo3.aov, 5)
plot(modelo3.aov, 6)
par(mfrow=c(1,1))
#### Box-Cox plot PHA
library(MASS)
bc <- boxcox(modelo3.aov)
(trans <- bc$x[which.max(bc$y)]) # conlambda < -2 o > 2 no debería usarse Box-Cox## [1] 0.5454545mnew <- lm(PHA^trans ~ SUELO * TRAT + REP:SUELO + REP:TRAT, data = hojasTRANS)
# QQ-plot PHA
op <- par(pty = "s", mfrow = c(1, 2))
qqnorm(modelo3.aov$residuals); qqline(modelo3.aov$residuals)
qqnorm(mnew$residuals); qqline(mnew$residuals)
par(op)
#### ANOVA para ARCSEN (opción 2)
modelo4.aov <- aov(ARCSEN ~ SUELO * TRAT + REP:SUELO + REP:TRAT, data = hojasTRANS) # transformados
summary(modelo3.aov)## Df Sum Sq Mean Sq F value Pr(>F)
## SUELO 1 582 582 6.124 0.0199 *
## TRAT 3 35783 11928 125.470 5.77e-16 ***
## SUELO:TRAT 3 186 62 0.653 0.5880
## SUELO:REP 2 201 101 1.059 0.3607
## TRAT:REP 3 214 71 0.750 0.5317
## Residuals 27 2567 95
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1shapiro.test(modelo4.aov$residuals)##
## Shapiro-Wilk normality test
##
## data: modelo4.aov$residuals
## W = 0.98233, p-value = 0.7749# Gráfico de diagnósticos
par(mfrow=c(1,2))
plot(modelo4.aov, 1)
plot(modelo4.aov, 2)
plot(modelo4.aov, 3)
plot(modelo4.aov, 4)
plot(modelo4.aov, 5)
plot(modelo4.aov, 6)
par(mfrow=c(1,1))
# Box-Cox plot ARCSEN
boxcox(modelo4.aov)
bc <- boxcox(modelo4.aov)
(trans <- bc$x[which.max(bc$y)])## [1] 0.7474747mnew <- lm(ARCSEN^trans ~ SUELO * TRAT + REP:SUELO + REP:TRAT, data = hojasTRANS)
# QQ-plot ARCSEN
op <- par(pty = "s", mfrow = c(1, 2))
qqnorm(modelo4.aov$residuals); qqline(modelo4.aov$residuals)
qqnorm(mnew$residuals); qqline(mnew$residuals)
par(op)
#### Comparación de pares de medias empleando la prueba Tukey
# PHA para SUELO
HSD.test.SUELO <- with(hojasTRANS, HSD.test(PHA, SUELO, DFerror = 13, MSerror = 46))
HSD.test.SUELO# ARSEN para SUELO
HSD.test.SUELO <- with(hojasTRANS, HSD.test(ARCSEN, SUELO, DFerror = 13, MSerror = 46))
HSD.test.SUELO# PHA para TRAT
HSD.test.TRAT <- with(hojasTRANS, HSD.test(PHA, TRAT, DFerror = 13, MSerror = 46))
HSD.test.TRAT# ARCSEN para TRAT
HSD.test.TRAT <- with(hojasTRANS, HSD.test(ARCSEN, TRAT, DFerror = 13, MSerror = 46))
HSD.test.TRAT
Códigos
en SAS
/***********************************************************
EXPERIMENTACIÓN EN AGRICULTURA
Fernández
Escobar, R.; Trapero, A.; Domínguez, J. 2010
CAPÍTULO 12:
DISEÑO EN BLOQUES DIVIDIDOS
************************************************************/
DATA
hojas;
INPUT
SUELO TRAT REP PHA;
LABEL SUELO = 'Sistema de
manejo del suelo'
TRAT =
'Momentos de los
tratamientos'
REP =
'Bloques'
PHA =
'Hojas enfermas (%)';
CARDS;
1 1 1 85.5
1 1 2 87.8
1 1 3 81.4
1 1 4 74.3
1 1 5 94.5
2 1 1 89.3
2 1 2 76.8
2 1 3 70.0
2 1 4 78.5
2 1 5 81.4
1 2 1 16.7
1 2 2 12.0
1 2 3 38.1
1 2 4 22.2
1 2 5 18.0
2 2 1 6.90
2 2 2 24.8
2 2 3 10.0
2 2 4 14.1
2 2 5 9.20
1 3 1 65.5
1 3 2 60.7
1 3 3 38.2
1 3 4 78.7
1 3 5 56.9
2 3 1 52.3
2 3 2 68.7
2 3 3 41.5
2 3 4 28.1
2 3 5 38.6
1 4 1 6.20
1 4 2 7.70
1 4 3 5.40
1 4 4 3.60
1 4 5 15.0
2 4 1 0.70
2 4 2 6.30
2 4 3 11.2
2 4 4 2.90
2 4 5 4.50
;
ODS HTML;
* datos originales;
PROC PRINT DATA =
hojas;
RUN;
* etiquetar variables categóricas;
PROC FORMAT;
VALUE SUELOCODIGO
1 = 'Cubierta
Vegetal'
2 = 'Suelo
Desnudo';
VALUE TRATCODIGO
1 = 'Sin
tratar'
2 = 'Primavera'
3 = 'Otoño'
4 = 'Prim
+ Otoño';
RUN;
DATA hojasID;
SET hojas;
FORMAT SUELO SUELOCODIGO. TRAT TRATCODIGO. REP
PHA;
RUN;
* datos con ID ;
PROC PRINT DATA =
hojasID;
RUN;
* transformando datos y agregando variable;
DATA hojasTRANS;
SET hojasID; *copiando datos de 'hojas';
ARCSEN =
ARSIN(SQRT(PHA/100))*180/constant("pi"); *transformar 'PHA' a ARSCEN;
LABEL = '% a Grados';
RUN;
* datos con variable transformada;
PROC PRINT DATA =
hojasTRANS;
RUN;
* estadísticas descriptivas;
PROC MEANS DATA =
hojasTRANS MEAN STD STDERR; *interacción
SUELO*TRAT;
VAR PHA
ARCSEN; CLASS SUELO TRAT; RUN;
PROC MEANS DATA =
hojasTRANS MEAN STD STDERR; *SUELO;
VAR
PHA; CLASS SUELO; RUN;
PROC MEANS DATA =
hojasTRANS MEAN STD STDERR; *TRAT;
VAR
PHA; CLASS TRAT; RUN;
PROC MEANS DATA =
hojasTRANS MEAN STD STDERR; *SUELO;
VAR
ARCSEN; CLASS SUELO; RUN;
PROC MEANS DATA =
hojasTRANS MEAN STD STDERR; *TRAT;
VAR
ARCSEN; CLASS TRAT; RUN;
* Prueba de Bartlett;
PROC ANOVA DATA = hojasTRANS;
CLASS
TRAT;
MODEL PHA
= TRAT;
MEANS
TRAT / HOVTEST = BARTLETT;
RUN;
PROC ANOVA DATA = hojasTRANS;
CLASS
TRAT;
MODEL
ARCSEN = TRAT;
MEANS
TRAT / HOVTEST = BARTLETT;
RUN;
* Opción 1: ANOVA bloques divididos;
PROC GLM DATA =
hojasTRANS;
CLASS SUELO TRAT REP;
MODEL ARCSEN = REP SUELO REP*SUELO TRAT REP*TRAT
SUELO*TRAT / SS2; *error tipo 2;
TEST H=REP SUELO E = REP*SUELO;
TEST H=REP
TRAT E =
REP*TRAT;
LSMEANS SUELO / STDERR E =
REP*SUELO;
LSMEANS
TRAT / STDERR E =
REP*TRAT;
LSMEANS
SUELO*TRAT/ ETYPE=2;
RUN;
* Opción 2: ANOVA
bloques divididos;
PROC MIXED DATA =
hojasTRANS METHOD = TYPE2;
CLASS
SUELO TRAT REP;
MODEL
ARCSEN = REP SUELO TRAT SUELO*TRAT;
RANDOM
REP*SUELO REP*TRAT;
LSMEANS
SUELO TRAT / PDIFF;
RUN;
* Comparación de medias PHA para TRAT;
PROC GLM;
CLASS REP
TRAT;
MODEL PHA
= REP TRAT REP*TRAT;
MEANS TRAT / TUKEY;
RUN;
* Comparación de medias ARCSEN para TRAT;
PROC GLM;
CLASS REP
TRAT;
MODEL
ARCSEN = REP TRAT REP*TRAT;
MEANS
TRAT / TUKEY;
RUN;
QUIT;
ODS HTML CLOSE;
----------------------------------------------
Si
desean adaptar los códigos de R para otros análisis, bajen tres archivos: tabla 12.1, capítulo 12.R, summarySE.R; para el análisis en SAS, el capítulo 12.sas. Finalmente, para mostrar cómo transformar una variable con arco-seno en Excel, incluimos tabla 12.1 – 12.2. El resumen y las referencias citadas en esta publicación proveen una discusión sobre la transformación.
No hay comentarios:
Publicar un comentario