#++++++++++++++++++++++++++++#
# UNIVARIATE APPROACH #
#++++++++++++++++++++++++++++#
#!____________________________
#! Qualitative Factor(s) (QL)
#!____________________________
#! Completely Randomized Design (CRD)
#! 1 factor - CRD - QLF
# Nonsense(experimental error = 0)
# Yi = mu + fe + e
r <- 2 # (repet. number)
fln <- 3 # (factor levels number)
crd00 <- gexp(mu = 0,
r = r,
fe = list(f1 = c(1, 2, 3)),
err = matrix(0,
nrow = r*fln),
round = 0)
crd00$X
print(crd00)
summary(crd00)
str(crd00)
#! 1 factor - CRD - QL
# Nonsense(error is 0)
# Yi = mu + fe + e
r <- 3 # (repet. number)
fln <- 5 # (factor levels number)
crd01 <- gexp(mu = 1,
r = r,
fe = list(f1 = c(0, 2, 4, 6, 8)),
err = matrix(0,
nrow = r*fln),
round = 2)
summary(crd01)
#! 1 factor - CRD - QL
# Default error: rmvnorm(sigma = diag(ncol(as.matrix([[fe]]))))
crd_1f <- gexp(mu = 1,
r = 3,
fe = list(f1 = c(1, 1, 5, 1, 1)),
fl = list(Treat = LETTERS[1:5]),
round = 2)
crd_1f$X
summary(crd_1f)
#! Binomial error - CRD - QL
e_binom <- as.matrix(rbinom(n = 15,
size = 5,
prob = 0.1))
crd_bin <- gexp(mu = 20,
err = e_binom,
r = 5,
fe = list(f1 = c(1, 4, 1)))
summary(crd_bin)
mod <- aov(Y1 ~ X1,
data = crd_bin$dfm)
shapiro.test(mod$res)
#! Factorial Experiment (FE) - CRD - QL
fe_crd00 <- gexp(mu = 0,
r = 2,
fe = list(f1 = c(1, 1, 5),
f2 = c(1, 1),
f3 = c(2, 2, 1)),
fl = list(A = paste('a',
1:3,
sep = ''),
B = paste('b',
1:2,
sep = ''),
C = paste('c',
1:3,
sep = '')),
round = 0,
type = 'FE')
fe_crd00$X
summary(fe_crd00)
#! Factorial Experiment (FE) - With interaction - CRD - QL
fe_crd01 <- gexp(mu = 30,
fe = list(f1 = c(1, 1, 3),
f2 = c(1, 1)),
fl = list(A = paste('a',
1:3,
sep = ''),
B = paste('b',
1:2,
sep = '')),
inte = c(3, 1, 1, 1, 1, 5), # (3*2)
round = 1,
type = 'FE')
summary(fe_crd01)
#! Split-plot Experiment (SPE) - CRD - QL
split_crd <- gexp(mu = 30,
fe = list(f1 = c(1, 1),
f2 = c(2, 3)),
fl = list(P = paste('p',
1:2,
sep = ''),
SP = paste('sp',
1:2,
sep = '')),
inte = c(1, 15, 1, 1), # (2*2)
round = 1,
type = 'SPE',
design = 'CRD')
split_crd$X
split_crd$Z
summary(split_crd)
split_crd01 <- gexp(mu = 30,
r = 3,
fe = list(f1 = c(1, 1),
f2 = c(2, 3),
f3 = c(1, 1, 1)),
fl = list(P = paste('p',
1:2,
sep = ''),
A = paste('a',
1:2,
sep = ''),
B = paste('b',
1:3,
sep = '')),
round = 1,
type = 'SPE',
design = 'CRD')
split_crd01$X
split_crd01$Z
summary(split_crd01)
#! Randomized Complete Block Design (RCBD) - QL
# 1 factor, 3 blocks
rcbd <- gexp(mu = 0,
r = 2,
fe = list(f1 = c(5, 1, 1)),
fl = list(TR = LETTERS[1:3]),
blke = c(1, 2, 3),
blkl = list(BLK = paste('B',
1:3,
sep = '')),
round = 1,
design = 'RCBD')
rcbd$X
summary(rcbd)
#! Factorial Experiment (FE) - RCBD - QL
fe_rcbd <- gexp(mu = 30,
r = 2,
fe = list(f1 = c(1, 1, 1),
f2 = c(2, 3)),
blke = c(1, 3),
inte = c(1, 15, 1, 1, 5, 1), # (3*2)
round = 1,
type = 'FE',
design = 'RCBD')
summary(fe_rcbd)
#! Multivariated - RCBD - QL
rcbd_m <- gexp(mu = c(0, 2),
fe = list(f1 = matrix(c(1, 1,
5, 1,
1, 1),
ncol = 2,
byrow = TRUE)),
blke = matrix(c(2, 1,
1, 2,
1, 1),
ncol = 2,
byrow = TRUE),
round = 1,
design = 'RCBD')
summary(rcbd_m)
#! Split-plot Experiment (SPE) - RCBD - QL
split_rcbd <- gexp(mu = 30,
r = 2,
fe = list(f1 = c(1, 1),
f2 = c(2, 3),
f3 = c(1, 1, 1)),
fl = list(P = paste('p',
1:2,
sep = ''),
B = paste('b',
1:2,
sep = ''),
C = paste('c',
1:3,
sep = '')),
blke = c(1, 2),
blkl = list(BLK = paste('B',
1:2,
sep = '')),
round = 1,
type = 'SPE',
design = 'RCBD')
split_rcbd$Z
summary(split_rcbd)
#! Latin Square Design (LSD) - QL
#!. Warning!!!! r = 5 by default
lsd00 <- gexp(design = 'LSD')
#Set r = 1 to omiting warning
lsd01 <- gexp(mu = 30,
r = 1,
fe = list(f1 = c(1, 1, 10)),
rowe = c(1, 1, 1),
cole = c(1, 1, 1),
rowl = list(Row = paste('r',
1:3,
sep = '')),
coll = list(Col = paste('c',
1:3,
sep = '')),
round = 0,
design = 'LSD')
summary(lsd01)
#! Factorial Experiment (FE) - LSD - QL
fe_lsd <- gexp(mu = 30,
r = 1,
fe = list(f1 = c(1, 1),
f2 = c(2, 3)),
rowe = c(1, 3, 2, 1),
cole = c(2, 2, 1, 1),
rowl = list(Row = paste('r',
1:4,
sep = '')),
coll = list(Col = paste('c',
1:4,
sep = '')),
inte = c(1, 15, 1, 1), # (2*2)
round = 1,
type = 'FE',
design = 'LSD')
summary(fe_lsd)
#! Split-plot Experiment (SPE) - LSD - QL
split_lsd <- gexp(mu = 30,
r = 1,
fe = list(f1 = c(1, 1, 2),
f2 = c(2, 3, 1)),
fl = list(P = paste('p',
1:3,
sep = ''),
SP = paste('sp',
1:3,
sep = '')),
inte = c(1, 15, 1, 1, 1, 1, 1, 1, 1), # (3*3)
rowe = c(1, 1, 1),
cole = c(1, 1, 1),
rowl = list(Row = paste('r',
1:3,
sep = '')),
coll = list(Col = paste('c',
1:3,
sep = '')),
round = 1,
type = 'SPE',
design = 'LSD')
summary(split_lsd)
#!_____________________________
#! Quantitative Factor(s) (QT)
#!_____________________________
#! CRD - Orthogonal polynomials
# Linear effect
# Nonsense(error is 0)
# Default contrasts: Orthogonal contrasts
r <- 4 # (repet. number)
fln <- 5 # (factor levels number)
level <- c(0, 10, 20, 30, 40)
crd_lo <- gexp(mu = 1, #in this case, mu=beta0 (intercept)
r = r,
fe = list(f1 = c(2, 0, 0, 0)), #b1 #b2 #b3 #b4
fl = list(Dose = level),
err = matrix(0,
nrow = r*fln),
round = 2)
crd_lo$X
summary(crd_lo)
plot(Y1 ~ Dose,
crd_lo$dfm)
# Quadratic effect
crd_qo <- gexp(mu = 2,
r = r,
fe = list(f1 = c(0, 3, 0, 0)), #b1 #b2 #b3 #b4
fl = list(Dose = level),
err = matrix(0,
nrow = r*fln))
summary(crd_qo)
plot(Y1 ~ Dose,
crd_qo$dfm)
# Cubic effect
crd_co <- gexp(mu = 2,
r = r,
fe = list(f1 = c(1, 1, 3, 0)), #b1 #b2 #b3 #b4
fl = list(Dose = level),
err = matrix(0,
nrow = r*fln))
summary(crd_co)
plot(Y1 ~ Dose,
crd_co$dfm)
# Not orthogonal polynomials
# Linear
cont_crd <- matrix(c(level,
level^2,
level^3,
level^4),
ncol = 4)
crd_l <- gexp(mu = 2,
r = 2,
fe = list(f1 = c(10, 0, 0, 0)), #b1 #b2 #b3 #b4
fl = list(Dose = level),
contrasts = list(Dose = cont_crd))
crd_l$X
summary(crd_l)
plot(Y1 ~ Dose,
crd_l$dfm)
reg <- lm(Y1 ~ Dose + I(Dose^2) + I(Dose^3) + I(Dose^4),
data = crd_l$dfm)
summary(reg)
# Linear and quadratic
level1 <- seq(0,30,by = 10)
cont_crd1 <- matrix(c(level1,
level1^2,
level1^3),
ncol = 3)
level2 <- 1:4
cont_crd2 <- matrix(c(level2,
level2^2,
level2^3),
ncol = 3)
crd_lq <- gexp(mu = 1,
r = 2,
fe = list(f1 = c(10, 0, 0), #b1 #b2 #b3
f2 = c(1, 8, 0)),
fl = list(P = level1,
N = level2),
contrasts = list(N = cont_crd2,
P = cont_crd1))
crd_lq$X
summary(crd_lq)
with(crd_lq$dfm,
plot(Y1 ~ P))
with(crd_lq$dfm,
plot(Y1 ~ N))
# Multivariated
crd_m <- gexp(mu = c(2, 10),
r = 4,
fe = list(f1 = matrix(c(10, 0, #L Q
0, 10,
0, 0),
ncol = 2,
byrow = TRUE)),
fl = list(Dose = level1),
contrasts = list(Dose = cont_crd1))
with(crd_m$dfm,
plot(Y1 ~ Dose))
with(crd_m$dfm,
plot(Y2 ~ Dose))
# RCBD - Orthogonal polynomios
level3 <- c(0, 2, 4, 6)
rcbd <- gexp(mu = 1,
fe = list(f1 = c(3, 0, 0)), #b1 #b2 #b3
blke = c(1, 2, 3),
r = 2,
fl = list(Dose = level3),
blkl = list(Blk = c('B1', 'B2', 'B3')),
design = 'RCBD')
rcbd$X
summary(rcbd)
# Not orthogonal
cont_crd3 <- matrix(c(level3, level3^2, level3^3),
ncol = 3)
rcbd_01 <- gexp(mu = 1,
fe = list(f1 = c(3, 0, 0)), #b1 #b2 #b3
blke = c(1, 2, 3),
r = 2,
fl = list(Dose = level3),
blkl = list(Blk = c('B1', 'B2', 'B3')),
contrasts = list(Dose = cont_crd3),
design = 'RCBD')
rcbd_01$X
summary(rcbd_01)
# Orthogonal polynomios - LSD
lsd <- gexp(mu = 1,
r = 1,
fe = list(f1 = c(3, 0, 0)), #b1 #b2 #b3
rowe = rep(1, 4),
cole = rep(1, 4),
fl = list(Dose = level1),
design = 'LSD')
lsd$X
summary(lsd)
lsd_01 <- gexp(mu = 1,
r = 1,
fe = list(f1 = c(3, 0, 0)), #b1 #b2 #b3
rowe = rep(1, 4),
cole = rep(1, 4),
rowl = list(row = paste('r',
1:4,
sep = '')),
fl = list(Dose = level1),
design = 'LSD')
lsd_01$X
summary(lsd_01)
# Not orthogonal
lsd_02 <- gexp(mu = 1,
r = 1,
fe = list(f1 = c(3, 0, 0)), #b1 #b2 #b3
rowe = rep(1, 4),
cole = rep(1, 4),
fl = list(Dose = level3),
contrasts = list(Dose = cont_crd3),
design = 'LSD')
lsd_02$X
str(lsd_02)
#!__________________________________________________________________________
#! Hibrid: qualitative and quantitative factors in the same experiment - HB
#!__________________________________________________________________________
#! CRD - HB
r <- 2 # (repet. number)
fl1 <- 4# (first factor levels number)
fl2 <- 3# (second factor levels number)
crd_hb <- gexp(mu = 1, #in this case, mu=beta0 (intercept)
r = r,
fe = list(f1 = c(2, 0, 0), #b1 #b2 #b3
f2 = c(1, 1, 3)),
fl = list(Dose = seq(0,30,
by = 10),
Trat = LETTERS[1:3]),
err = matrix(0,
nrow = r*fl1*fl2),
round = 2)
crd_hb$X
summary(crd_hb)
#Only one contrasts!
crd_hb2 <- gexp(mu = 1, #in this case, mu=beta0 (intercept)
r = r,
fe = list(f1 = c(2, 0, 0), #b1 #b2 #b3
f2 = c(1, 1, 3)),
fl = list(Dose = level1,
Trat = LETTERS[1:3]),
err = matrix(0,
nrow = r*fl1*fl2),
contrasts = list(Dose = cont_crd1),
round = 2)
crd_hb2$X
summary(crd_hb)
#! RCBD - HB
r <- 2
blke <- c(1, 2)
level <- c(0, 10, 20, 30)
(error <- matrix(rep(0,
4^1*3^1*r*length(blke)),
ncol=1))
rcbd_hb <- gexp(mu = 2,
err = error,
r = r,
fe = list(f1 = c(0, 1, 0), # Qualitative
f2 = c(1, 0, 0)), # Quantitative linear
fl = list(Var = LETTERS[1:3],
Dose = level),
blke = blke,
blkl = list(Blk = c('B1', 'B2')),
design = 'RCBD')
rcbd_hb$X
summary(rcbd_hb)
str(rcbd_hb)
#! LSD - QT
set.seed(3)
lsd <- gexp(mu = 100,
r = 1,
fe = list(f1 = c(10, # b1
20, # b2
0, # b3
0)), # b4
fl = list(tra = seq(0,
40,
by = 10)),
rowe = c(1, 2, 3, 4, 5),
rowl = list(row = paste('r',
1:5,
sep = '')),
cole = c(5, 4, 3, 2, 1),
coll = list(col = paste('c',
1:5,
sep = '')),
design = 'LSD')
summary(lsd)
plot(Y1 ~ tra, lsd$dfm)
#! FE - LSD - QT
fe_lsd <- gexp(mu = 10,
fe = list(f1 = c(2, 3),
f2 = c(5, # b1*
0, # b2
0, # b3
0)), # b4
rowe = rep(1, 10),
cole = rep(1, 10),
fl = list(var = paste('v',
1:2,
sep = ''),
tra = seq(0,
40,
by = 10)),
coll = list(col = paste('c',
1:10,
sep = '')),
rowl = list(row = paste('r',
1:10,
sep = '')),
type = 'FE',
design = 'LSD')
fe_lsd$X
str(fe_lsd)
summary(fe_lsd)
plot(Y1 ~ tra,
fe_lsd$dfm)
#! SPE - QL - QT
spe_lsd <- gexp(mu = 100,
r = 1,
fe = list(f1 = c(2, 3, 1),
f2 = c(1, # b1
5, # b2*
1)), # b3
fl = list(p = paste('p',
1:3,
sep = ''),
sp = seq(0,
30,
by = 10)),
rowe = c(1, 2, 3),
cole = c(3, 2, 1),
rowl = list(row = paste('r',
1:3,
sep = '')),
coll = list(col = paste('c',
1:3,
sep = '')),
round = 1,
type = 'SPE',
design = 'LSD')
summary(spe_lsd)
plot(spe_lsd)
#++++++++++++++++++++++++++++#
# MULTIVARIATE APPROACH #
#++++++++++++++++++++++++++++#
#! CRD - QL
# Error = 0 - Nonsense (you can easily undertand the effects)
r <- 2 # (repet. number)
fln <- 3 # (factor levels number)
crd_m01 <- gexp(mu = c(0,10),
r = r,
fe = list(f1 = matrix(c(1, 0, #Y1 Y2
2, 1,
3, 3),
ncol = 2,
byrow = TRUE)),
err = mvtnorm::rmvnorm(n = fln * r,
sigma = matrix(c(0, 0,
0, 0),
ncol = 2)),
round = 0)
summary(crd_m01)
#! FE - CRD - QL
r <- 2
crd_m_fe01 <- gexp(mu = c(0, 0),
r = r,
err = mvtnorm::rmvnorm(n = 3^1 * 2^1 * r,
sigma = matrix(c(0, 0,
0, 0),
ncol = 2)),
fe = list(f1 = matrix(c(0, 3, #X1 X1
1, 4, #X2 X2
2, 5), #X3 X3
ncol = 2,
byrow = TRUE),
f2 = matrix(c(0, 2, #X1 X1
1, 3), #X2 X2
ncol = 2,
byrow = TRUE)),
type = 'FE',
round = 1)
summary(crd_m_fe01)
#! FE - CRD - QL
# Using default error
set.seed(30)
crd_m_fe02 <- gexp(mu = c(0, 2),
r = 3,
fe = list(f1 = matrix(c(1, 1,
5, 1,
1, 1),
ncol = 2,
byrow = TRUE),
f2 = matrix(c(1, 3,
2, 2),
ncol = 2,
byrow = TRUE)),
type = 'FE',
round = 1)
summary(crd_m_fe02)
#! SPE - CRD - QL
# Using default error
crd_m_spe01 <- gexp(mu = c(0, 2),
r = 3,
fe = list(f1 = matrix(c(1, 1,
5, 1,
1, 1),
ncol = 2,
byrow = TRUE),
f2 = matrix(c(1, 3,
2, 2),
ncol = 2,
byrow = TRUE)),
type = 'SPE',
round = 1)
summary(crd_m_spe01)
#! RCBD - QL
r <- 2 # (repet. number)
fln <- 3 # (factor levels number)
bln <- 3 # (block levels number)
rcbd_m01 <- gexp(mu = c(0,10),
r = r,
fe = list(f1 = matrix(c(1, 0, #Y1 Y2
2, 1,
3, 3),
ncol = 2,
byrow = TRUE)),
blke = matrix(c(2, 1,
4, 1,
6, 1),
ncol = 2,
byrow = TRUE),
err = mvtnorm::rmvnorm(n = fln * r * bln,
sigma = matrix(c(0, 0,
0, 0),
ncol = 2)),
design = 'RCBD',
round = 0)
summary(rcbd_m01)
#! FE - RCBD - QL
rcbd_m_fe01 <- gexp(mu = c(0, 0),
r = r,
err = mvtnorm::rmvnorm(n = 3^1 * 2^1 * r * bln,
sigma = matrix(c(0, 0,
0, 0),
ncol = 2)),
fe = list(f1 = matrix(c(0, 3, #X1 X1
1, 4, #X2 X2
2, 5), #X3 X3
ncol = 2,
byrow = TRUE),
f2 = matrix(c(0, 2, #X1 X1
1, 3), #X2 X2
ncol = 2,
byrow = TRUE)),
blke = matrix(c(2, 1,
4, 1,
6, 1),
ncol = 2,
byrow = TRUE),
type = 'FE',
design = 'RCBD',
round = 1)
summary(rcbd_m_fe01)
#! SPE - RCBD - QL
rcbd_m_spe01 <- gexp(mu = c(0, 2),
r = 2,
fe = list(f1 = matrix(c(1, 1,
5, 1,
1, 1),
ncol = 2,
byrow = TRUE),
f2 = matrix(c(1, 3,
2, 2),
ncol = 2,
byrow = TRUE),
f3 = matrix(c(1, 3,
2, 2),
ncol = 2,
byrow = TRUE)),
blke = matrix(c(2, 1,
4, 1,
6, 1),
ncol = 2,
byrow = TRUE),
type = 'SPE',
design = 'RCBD',
round = 1)
summary(rcbd_m_spe01)
#! LSD - QL
lsd_m01 <- gexp(mu = c(0,10),
r = 1,
fe = list(f1 = matrix(c(1, 0,
2, 1,
3, 3),
ncol = 2,
byrow = TRUE)),
rowe = matrix(rep(1, 6),
ncol = 2),
cole = matrix(rep(1, 6),
ncol = 2),
err = mvtnorm::rmvnorm(n = 3^2,
sigma = matrix(c(0, 0,
0, 0),
ncol = 2)),
design = 'LSD',
round = 0)
summary(lsd_m01)
#! LSD/FE - QL
lsd_m_fe01 <- gexp(mu = c(0, 0),
r = 1,
err = mvtnorm::rmvnorm(n = 3^1 * 2^1 * 6,
sigma = matrix(c(0, 0,
0, 0),
ncol = 2)),
#Y1 Y2
fe = list(f1 = matrix(c(0, 3, #X1 X1
1, 4, #X2 X2
2, 5), #X3 X3
ncol = 2,
byrow = TRUE),
#Y1 Y2
f2 = matrix(c(0, 2, #X1 X1
1, 3), #X2 X2
ncol = 2,
byrow = TRUE)),
rowe = matrix(rep(1, 12),
ncol = 2),
cole = matrix(rep(1, 12),
ncol = 2),
type = 'FE',
design = 'LSD',
round = 1)
summary(lsd_m_fe01)
#! SPE - LSD - QL
# Using default error
lsd_m_spe01 <- gexp(mu = c(0, 2),
r = 1,
fe = list(f1 = matrix(c(1, 1,
5, 1,
1, 1),
ncol = 2,
byrow = TRUE),
f2 = matrix(c(1, 3,
2, 2),
ncol = 2,
byrow = TRUE)),
rowe = matrix(rep(1, 6),
ncol = 2),
cole = matrix(rep(1, 6),
ncol = 2),
type = 'SPE',
design = 'LSD',
round = 1)
summary(lsd_m_spe01)
#! FE - RCBD - QL
r = 1
bln = 3
fe_rcbd_m <- gexp(mu = c(0, 0),
r = 1,
err = mvtnorm::rmvnorm(n = 3^1 * 2^1 * r * bln,
sigma = matrix(c(0, 0,
0, 0),
ncol = 2)),
fe = list(f1 = matrix(c(0, 3, #X1 X1
1, 4, #X2 X2
2, 5), #X3 X3
ncol = 2,
byrow = TRUE),
f2 = matrix(c(0, 2, #X1 X1
1, 3), #X2 X2
ncol = 2,
byrow = TRUE)),
blke = matrix(c(2, 1,
4, 1,
6, 1),
ncol = 2,
byrow = TRUE),
type = 'FE',
design = 'RCBD')
str(fe_rcbd_m)
summary(fe_rcbd_m)
#! SPE - RCBD - QL
spe_rcbd_m <- gexp(mu = c(0, 2),
r = 3,
fe = list(f1 = matrix(c(1, 1,
5, 1,
1, 1),
ncol = 2,
byrow = TRUE),
f2 = matrix(c(1, 3,
2, 2),
ncol = 2,
byrow = TRUE),
f3 = matrix(c(1, 3,
2, 2),
ncol = 2,
byrow = TRUE)),
blke = matrix(c(2, 1,
4, 1,
6, 1),
ncol = 2,
byrow = TRUE),
type = 'SPE',
design = 'RCBD')
str(spe_rcbd_m)
summary(spe_rcbd_m)
Run the code above in your browser using DataLab