CE 074 - Controle de Processos Industriais

http://www.leg.ufpr.br/ce074


Universidade Federal do Paraná Prof. Fernando de Pol Mayer

Curso de Graduação em Estatística Lab. de Estatística e Geoinformação - LEG

Departamento de Estatística - UFPR


Técnicas de confundimento para blocagem em experimentos fatoriais \(2^k\)

Existem certas situações em que é praticamente impossível fazer todas as corridas de um experimento em condições uniformes. Por exemplo, pode haver limitações da quatidade de matéria prima, ou matéria prima de diversas origens. As condições de contorno podem mudar ao longo do ensaio (temperatura, ventilação, luz). Pode existir um número elevado de tratamentos difícil de acomodar em um curto espaço de tempo/espaço ou reduzído número de instrumentos/operadores, além de ser desejável variar as condições de contorno para garantir eficiência/robustez aos resultados. A técnica experimental adotada nessas situações é a blocagem.

A idéia central da blocagem é fazer com que as unidades experimentais (UEs) sejam homogêneas dentro dos blocos. Os blocos são completos quando em cada bloco existe pelo menos uma UE de cada tratamento, e incompleto caso contrário.

Nos experimentos \(2^k\) existe uma série de opções de blocagem. A primeira é repetir o experimento de forma que cada repetição completa (que inclui todos os tratamentos) seja um bloco. É o caso comum quando tem-se poucos tratamentos (geralmente \(2^2\) ou \(2^3\)), e nesses casos específicos temos um fatorial com blocos completos.

Como nos experimentos fatorias \(2^k\) (\(k \geq 3\)) o número de tratamentos geralmente é grande, devido ao caráter exploratório do experimento, os blocos dificilmente cumprirão seu papel se forem completos, por isso geralmente adota-se blocos incompletos. Nesse caso os tratamentos devem ser particionados e atribuídos aos blocos. Nada impede que essa partição dos tratamentos seja aleatória, porém quando feita estrategicamente leva-se algumas vantagens.

A estratégia adotada para se atribuir os tratamentos aos blocos é a de confundimento. A idéia central é tomar interações de alta ordem e propositalmente confundir o efeito dessa interação com o efeito dos blocos. Isso porque interações de ordem alta dificilmente são interpretáveis, e o efeito dos blocos não é do interesse do pesquisador. O bloco está presente para acomodar alterações das condições de contorno. Dessa forma, não é desconforto ter esses efeitos confundidos/misturados quando o foco do experimento são os efeitos principais e interações de ordem mais baixa.

Vamos considerar a construção e análise de fatoriais \(2^k\) em \(2^p\) blocos incompletos, onde \(p < k\). Consequentemente, estes experimentos podem ser divididos em 2, 4, 8, \(\ldots\) blocos.

Blocagem em um experimento fatorial \(2^k\) com repetição

Se um experimento fatorial \(2^k\) for replicado \(n\) vezes sob condições não homogêneas, então cada conjunto destas condições com todos os tratamentos definem um bloco. Portanto, teríamos \(n\) blocos completos. As corridas em cada bloco devem ser realizadas de forma aleatória.

A análise do experimento nesse caso é idêntica àquela de um experimento fatorial \(2^k\) sem blocos, com a exceção de que haverá também uma soma de quadrados para bloco, dada por

\[ SQ_{Bloco} = \frac{\sum_{i=1}^n B_i^2}{2^k} - \frac{y_{...}^2}{N} \]

onde \(B_i\) é o total de cada bloco \(i\), \(y_{...}^2\) é o quadrado da soma total de todas as observações, e \(N\) é o número total de corridas do experimento.

Exemplo

Experimento \(2^2\) em blocos completos

Um processo químico é investigado em relação à dois fatores e o experimento é conduzido em 3 blocos completos. Os resultados estão abaixo. Faça a análise estatística dos resultados.

bloco 1  b 2  b 3
(1) = 28   25   27
  a = 36   32   32
  b = 18   19   23
 ab = 31   30   29
da <- expand.grid(A = c(-1, 1),
                  B = c(-1, 1),
                  bloco = c("I", "II", "III"))
da
    A  B bloco
1  -1 -1     I
2   1 -1     I
3  -1  1     I
4   1  1     I
5  -1 -1    II
6   1 -1    II
7  -1  1    II
8   1  1    II
9  -1 -1   III
10  1 -1   III
11 -1  1   III
12  1  1   III
da$y <- c(28,36,16,31,25,32,19,30,27,32,23,29)
## Note que bloco deve ser fator
str(da)
'data.frame':   12 obs. of  4 variables:
 $ A    : num  -1 1 -1 1 -1 1 -1 1 -1 1 ...
 $ B    : num  -1 -1 1 1 -1 -1 1 1 -1 -1 ...
 $ bloco: Factor w/ 3 levels "I","II","III": 1 1 1 1 2 2 2 2 3 3 ...
 $ y    : num  28 36 16 31 25 32 19 30 27 32 ...
 - attr(*, "out.attrs")=List of 2
  ..$ dim     : Named int  2 2 3
  .. ..- attr(*, "names")= chr  "A" "B" "bloco"
  ..$ dimnames:List of 3
  .. ..$ A    : chr  "A=-1" "A= 1"
  .. ..$ B    : chr  "B=-1" "B= 1"
  .. ..$ bloco: chr  "bloco=I" "bloco=II" "bloco=III"
## Ajustando o modelo
m0 <- lm(y ~ bloco + A * B, data = da)
anova(m0)
Analysis of Variance Table

Response: y
          Df  Sum Sq Mean Sq F value    Pr(>F)    
bloco      2   4.167   2.083  0.3304 0.7309303    
A          1 225.333 225.333 35.7357 0.0009834 ***
B          1  85.333  85.333 13.5330 0.0103463 *  
A:B        1  12.000  12.000  1.9031 0.2169434    
Residuals  6  37.833   6.306                      
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Ajustando o modelo sem a interação, que não foi significativa
m1 <- update(m0, . ~ . - A:B)
anova(m1)
Analysis of Variance Table

Response: y
          Df  Sum Sq Mean Sq F value    Pr(>F)    
bloco      2   4.167   2.083  0.2926 0.7549907    
A          1 225.333 225.333 31.6522 0.0007941 ***
B          1  85.333  85.333 11.9866 0.0105166 *  
Residuals  7  49.833   7.119                      
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Experimento \(2^3\) em blocos completos

Uma planejamento fatorial \(2^3\) foi corrido em um processo químico. Os fatores do planejamento são A = tempo, B = concentração, C = pressão. Duas repetições foram feitas em horas distintas do dia idetificadas pelo nível de bloco. A variável resposta é o rendimento. Os dados estão disponíveis com os comandos a seguir. Faça a análise estatística dos resultados.

da <- expand.grid(A = c(-1, 1),
                  B = c(-1, 1),
                  C = c(-1, 1),
                  bloco = c("I", "II"))
da$y <- c(12,18,13,16,17,15,20,25,10,25,13,24,19,21,17,23)
da
    A  B  C bloco  y
1  -1 -1 -1     I 12
2   1 -1 -1     I 18
3  -1  1 -1     I 13
4   1  1 -1     I 16
5  -1 -1  1     I 17
6   1 -1  1     I 15
7  -1  1  1     I 20
8   1  1  1     I 25
9  -1 -1 -1    II 10
10  1 -1 -1    II 25
11 -1  1 -1    II 13
12  1  1 -1    II 24
13 -1 -1  1    II 19
14  1 -1  1    II 21
15 -1  1  1    II 17
16  1  1  1    II 23
str(da)
'data.frame':   16 obs. of  5 variables:
 $ A    : num  -1 1 -1 1 -1 1 -1 1 -1 1 ...
 $ B    : num  -1 -1 1 1 -1 -1 1 1 -1 -1 ...
 $ C    : num  -1 -1 -1 -1 1 1 1 1 -1 -1 ...
 $ bloco: Factor w/ 2 levels "I","II": 1 1 1 1 1 1 1 1 2 2 ...
 $ y    : num  12 18 13 16 17 15 20 25 10 25 ...
 - attr(*, "out.attrs")=List of 2
  ..$ dim     : Named int  2 2 2 2
  .. ..- attr(*, "names")= chr  "A" "B" "C" "bloco"
  ..$ dimnames:List of 4
  .. ..$ A    : chr  "A=-1" "A= 1"
  .. ..$ B    : chr  "B=-1" "B= 1"
  .. ..$ C    : chr  "C=-1" "C= 1"
  .. ..$ bloco: chr  "bloco=I" "bloco=II"
## Ajuste do modelo
m0 <- lm(y ~ bloco + A * B * C, data = da)
anova(m0)
Analysis of Variance Table

Response: y
          Df Sum Sq Mean Sq F value  Pr(>F)   
bloco      1  16.00  16.000  1.6232 0.24331   
A          1 132.25 132.250 13.4167 0.00804 **
B          1  12.25  12.250  1.2428 0.30175   
C          1  42.25  42.250  4.2862 0.07718 . 
A:B        1   1.00   1.000  0.1014 0.75939   
A:C        1  36.00  36.000  3.6522 0.09760 . 
B:C        1   9.00   9.000  0.9130 0.37113   
A:B:C      1  20.25  20.250  2.0543 0.19489   
Residuals  7  69.00   9.857                   
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Ajustando o modelo apenas com efeitos significativos
m1 <- update(m0, . ~ bloco + A * C)
anova(m1)
Analysis of Variance Table

Response: y
          Df Sum Sq Mean Sq F value   Pr(>F)   
bloco      1  16.00  16.000  1.5785 0.235007   
A          1 132.25 132.250 13.0471 0.004083 **
C          1  42.25  42.250  4.1682 0.065921 . 
A:C        1  36.00  36.000  3.5516 0.086171 . 
Residuals 11 111.50  10.136                    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Confundimento do fatorial \(2^k\) em dois blocos

A idéia central é dividir as \(2^k\) UEs igualmente em dois blocos de forma que uma interação de ordem que não tenha interesse seja confundida com o efeito dos blocos. É usual usar a interação de maior ordem para fazer a divisão, que é simples: aquelas UEs com sinal (-) serão de um bloco e as de sinal (+) serão do outro. Pode-se fazer essa separação obtendo a coluna de sinais da interação mais alta ou usando o contraste de definição.

Blocagem de um fatorial \(2^2\) em dois blocos

Considere um planejamento \(2^2\) onde cada uma das 4 combinações de tratamentos requeira quatro horas de análise de laboratório. Dessa forma, dois dias são necessários para realizar o experimento. Se dias forem considerados como blocos, então temos que atribuir duas das quatro combinações em cada dia.

Este experimento está representado na figura abaixo

                                                     Bloco I        Bloco II
[+] b--------------ab                +
    |               |                |              +-------+      +-------+
    |               |                |              |       |      |       |
    |               |                |              |  (1)  |      |   a   |
 B  |               |                |              |       |      |       |
    |               |                |              |       |      |       |
    |               |                |              |  ab   |      |   b   |
    |               |                |              |       |      |       |
[-] (1)-------------a                +              +-------+      +-------+

   [-]      A      [+]

Note que o bloco I contém as combinações de tratamento (1) e ab, e que o bloco II contém a e b. Lembrando que os contrastes para estimar os efeitos dos fatores A e B são

\[ \begin{align} contr_A = ab + a - b - (1) \\ contr_B = ab + b - a - (1) \end{align} \]

Observe que estes contrastes não são afetados pela blocagem, uma vez que em cada contraste há uma combinação de tratamentos mais e outra menos, provenientes de cada bloco. Portanto, qualquer diferença entre o bloco I e o bloco II será cancelada.

O contraste para a interação é

\[ contr_{AB} = ab + (1) - a - b \]

Já que as duas combinações de tratamento com sinal mais (ab e (1)) estão no bloco I, e as duas com sinal menos estão no bloco II (a e b), o efeito do bloco e da interação AB é o mesmo. Ou seja, a interação AB está confundida com os blocos.

A razão para isso está clara quando analisamos a tabela de sinais para o planejamento \(2^2\)

     A  B AB
(1) -1 -1  1
a    1 -1 -1
b   -1  1 -1
ab   1  1  1

Dessa tabela, vemos que todas as combinações de tratamentos que possuem sinal mais em AB são atribuídas ao bloco I, enquanto que os tratamentos com sinal menos em AB são atribuídas ao bloco II.

Essa abordagem pode ser usada para confundir qualquer efeito (A, B, ou AB) com blocos. Por exemplo, se a e ab estivessem no bloco I, e (1) e b no bloco II, então seria o efeito A que estaria confundido com blocos. A prática usual é confundir a interação de ordem mais alta com blocos, já que geralmente essa interação não tem interpretação prática e normalmente também não é significativa.

Além disso, a definição de dois blocos em qualquer esquema fatorial \(2^k\) pode ser feita por essa abordagem.

Blocagem de um fatorial \(2^3\) em dois blocos

Considere um planejamento fatorial \(2^3\). Para definir os tratamentos que serão atribuídos a cada bloco, vamos considerar a interação de ordem mais alta, ABC. Pela tabela de sinais desse planejamento, vamos atribuir os tratamentos com sinal menos na coluna ABC ao bloco I, e os tratamentos com sinal mais na coluna ABC ao bloco II.

     A  B  C ABC
(1) -1 -1 -1  -1
a    1 -1 -1   1
b   -1  1 -1   1
ab   1  1 -1  -1
c   -1 -1  1   1
ac   1 -1  1  -1
bc  -1  1  1  -1
abc  1  1  1   1

O planejamento resultante pode ser visto na representação geométrica abaixo.

           Bloco I           |             Bloco II
                             |
        bc --------          |            -------- abc
        .|        .|         |          .|        .|
       __|______ab |         |        b__|______   |
      |  |      |  |         |        |  |      |  |
    B |   ------|--ac        |      B |  c -----|--
      | .       | . C        |        | .       | . C
     (1)________|            |        |_________a
           A                 |             A

Novamente, é importante lembrar que a combinação de tratamentos dentro de cada bloco deve ser atribuída de forma aleatória.

Blocagem de experimentos fatoriais \(2^k\) em dois blocos usando contraste de definição

Outro método mais geral para construir os blocos é através dos contrastes de definição. Este método usa uma combinção linear

\[ L = \alpha_1 x_1 + \alpha_2 x_2 + \cdots + \alpha_k x_k \]

onde \(x_i\) é o nível do \(i\)-ésimo fator aparecendo em uma combinação de tratamentos, e \(\alpha_i\) é o expoente aparecendo no \(i\)-ésimo fator no efeito que deve ser confundido. Por exemplo, se o efeito a ser confundido for ABCD, então os valores de \(\alpha_i\) serão todos iguais a 1, pois \(ABCD = A^1B^1C^1D^1\). Se o efeito a ser confundido for ACD, então os valores serão \(\alpha_1 = \alpha_3 = \alpha_4 = 1\), e \(\alpha_2 = 0\), pois \(ACD = A^1B^0C^1D^1\).

Dessa forma, para o sistema \(2^k\), temos tanto \(\alpha_i = 0\) ou \(1\), e \(x_i = 0\) (nível baixo) ou \(1\) (nível alto). Combinações de tratamentos que produzam o mesmo valor de \(L \mod 2\) serão colocados no mesmo bloco. Como os únicos valores possíveis de \(L \mod 2\) são 0 e 1, isso atribuirá as \(2^k\) combinações de tratamentos à extamente dois blocos.

NOTA: a função \(x \mod 2\) retorna o resto da divisão de x pelo número 2. \(\text{mod}\) é a função módulo, e no R é reprsentada por %%.

Como exemplo, considere um planejamento \(2^3\), com a interação ABC (a de ordem mais alta) confundida com bloco. Aqui, \(x_1\) corresponde a A, \(x_2\) a B, e \(x_3\) a C. Portanto, temos que \(\alpha_1 = \alpha_2 = \alpha_3 = 1\), pois como o fator a ser confundido é ABC, então o expoente destes três fatores é 1. Portanto, o contraste de definição utilizado para confundir ABC com blocos é

\[ L = x_1 + x_2 + x_3 \]

Com a finalidade de atribuir as combinações de tratamentos aos dois blocos, substituímos as combinações de tratamentos ao contraste de definição, como segue:

\[ \begin{align} (1):& \quad L = 1(0) + 1(0) + 1(0) = 0 \mod 2 = 0 \\ a:& \quad L = 1(1) + 1(0) + 1(0) = 1 \mod 2 = 1 \\ b:& \quad L = 1(0) + 1(1) + 1(0) = 1 \mod 2 = 1 \\ ab:& \quad L = 1(1) + 1(1) + 1(0) = 2 \mod 2 = 0 \\ c:& \quad L = 1(0) + 1(0) + 1(1) = 1 \mod 2 = 1 \\ ac:& \quad L = 1(1) + 1(0) + 1(1) = 2 \mod 2 = 0 \\ bc:& \quad L = 1(0) + 1(1) + 1(1) = 2 \mod 2 = 0 \\ abc:& \quad L = 1(1) + 1(1) + 1(1) = 3 \mod 2 = 1 \end{align} \]

NOTE que na notação \((0,1)\), a combinação (1) é representada por 000, a por 100, e assim por diante.

Dessa forma, as combinações (1), ab, ac, e bc são corridas no bloco I, enquanto que a, b, c, e abc são corridas no bloco 2. Veja que esta atribuição é idêntica àquela realizada ao se utilizar a coluna ABC da tabela de sinais. O contraste de definição é apenas uma generalização daquele método.

Exemplos

##======================================================================
## Fatorial 2^3 com a notação (0,1)
da <- expand.grid(A = c(0, 1),
                  B = c(0, 1),
                  C = c(0, 1))
row.names(da) <- c("(1)", "a", "b", "ab", "c", "ac", "bc", "abc")
da
    A B C
(1) 0 0 0
a   1 0 0
b   0 1 0
ab  1 1 0
c   0 0 1
ac  1 0 1
bc  0 1 1
abc 1 1 1
## Montando o contraste de definição, considerando a interação ABC
## confundida com bloco
## alpha = 1,1,1 pois A, B e C estão presentes nesse efeito
alpha <- c(1, 1, 1) 

## Cálculo dos contrastes para cada combinação de tratamento, já
## considerando mod 2
## (1)
sum(alpha * da[1, ]) %% 2
[1] 0
## a
sum(alpha * da[2, ]) %% 2
[1] 1
## E assim por diante... Para facilitar podemos usar apply
apply(da, 1, function(x) sum(alpha * x) %% 2)
(1)   a   b  ab   c  ac  bc abc 
  0   1   1   0   1   0   0   1 
## E criar uma função para calcular L mod 2
contr.def <- function(alpha, x) sum(alpha * x) %% 2
da$bloco <- apply(da, 1, contr.def, alpha = alpha)
da
    A B C bloco
(1) 0 0 0     0
a   1 0 0     1
b   0 1 0     1
ab  1 1 0     0
c   0 0 1     1
ac  1 0 1     0
bc  0 1 1     0
abc 1 1 1     1
##======================================================================
## Dividindo um fatorial 2^4 em dois blocos, pela tabela de sinais
da <- do.call(expand.grid, replicate(4, list(c(-1, 1))))
names(da) <- LETTERS[1:ncol(da)]
row.names(da) <- apply(da, 1,
                       function(i) paste(letters[1:4][i==1], collapse = ""))
row.names(da)[1] <- "(1)"
da
      A  B  C  D
(1)  -1 -1 -1 -1
a     1 -1 -1 -1
b    -1  1 -1 -1
ab    1  1 -1 -1
c    -1 -1  1 -1
ac    1 -1  1 -1
bc   -1  1  1 -1
abc   1  1  1 -1
d    -1 -1 -1  1
ad    1 -1 -1  1
bd   -1  1 -1  1
abd   1  1 -1  1
cd   -1 -1  1  1
acd   1 -1  1  1
bcd  -1  1  1  1
abcd  1  1  1  1
## Usando a interação de ordem mais alta, ABCD, para confundir com blocos
da$ABCD <- with(da, A * B * C * D)
da
      A  B  C  D ABCD
(1)  -1 -1 -1 -1    1
a     1 -1 -1 -1   -1
b    -1  1 -1 -1   -1
ab    1  1 -1 -1    1
c    -1 -1  1 -1   -1
ac    1 -1  1 -1    1
bc   -1  1  1 -1    1
abc   1  1  1 -1   -1
d    -1 -1 -1  1   -1
ad    1 -1 -1  1    1
bd   -1  1 -1  1    1
abd   1  1 -1  1   -1
cd   -1 -1  1  1    1
acd   1 -1  1  1   -1
bcd  -1  1  1  1   -1
abcd  1  1  1  1    1
da <- da[order(da$ABCD),]
da
      A  B  C  D ABCD
a     1 -1 -1 -1   -1
b    -1  1 -1 -1   -1
c    -1 -1  1 -1   -1
abc   1  1  1 -1   -1
d    -1 -1 -1  1   -1
abd   1  1 -1  1   -1
acd   1 -1  1  1   -1
bcd  -1  1  1  1   -1
(1)  -1 -1 -1 -1    1
ab    1  1 -1 -1    1
ac    1 -1  1 -1    1
bc   -1  1  1 -1    1
ad    1 -1 -1  1    1
bd   -1  1 -1  1    1
cd   -1 -1  1  1    1
abcd  1  1  1  1    1
##======================================================================
## Dividindo um fatorial 2^4 em dois blocos, usando o contraste de
## definição
db <- do.call(expand.grid, replicate(4, list(c(0, 1))))
names(db) <- LETTERS[1:ncol(db)]
row.names(db) <- apply(db, 1,
                       function(i) paste(letters[1:4][i==1], collapse = ""))
row.names(db)[1] <- "(1)"
db
     A B C D
(1)  0 0 0 0
a    1 0 0 0
b    0 1 0 0
ab   1 1 0 0
c    0 0 1 0
ac   1 0 1 0
bc   0 1 1 0
abc  1 1 1 0
d    0 0 0 1
ad   1 0 0 1
bd   0 1 0 1
abd  1 1 0 1
cd   0 0 1 1
acd  1 0 1 1
bcd  0 1 1 1
abcd 1 1 1 1
## Usando a interação de ordem mais alta, ABCD, para confundir com
## blocos. Dessa forma temos
## L = x_1 + x_2 + x_3 + x_4
## com alpha_i = 1
alpha <- c(1, 1, 1, 1)
db$bloco <- apply(db, 1, contr.def, alpha = alpha)
db
     A B C D bloco
(1)  0 0 0 0     0
a    1 0 0 0     1
b    0 1 0 0     1
ab   1 1 0 0     0
c    0 0 1 0     1
ac   1 0 1 0     0
bc   0 1 1 0     0
abc  1 1 1 0     1
d    0 0 0 1     1
ad   1 0 0 1     0
bd   0 1 0 1     0
abd  1 1 0 1     1
cd   0 0 1 1     0
acd  1 0 1 1     1
bcd  0 1 1 1     1
abcd 1 1 1 1     0
db[order(db$bloco), ]
     A B C D bloco
(1)  0 0 0 0     0
ab   1 1 0 0     0
ac   1 0 1 0     0
bc   0 1 1 0     0
ad   1 0 0 1     0
bd   0 1 0 1     0
cd   0 0 1 1     0
abcd 1 1 1 1     0
a    1 0 0 0     1
b    0 1 0 0     1
c    0 0 1 0     1
abc  1 1 1 0     1
d    0 0 0 1     1
abd  1 1 0 1     1
acd  1 0 1 1     1
bcd  0 1 1 1     1
## Para variar, vamos escolhar e interação tripla ACD para definir o
## contraste de forma a separar as corridas em 2 blocos. Nossa função de
## definição fica
## L = x_1 + x_3 + x_4
## Note que todos os coeficientes alpha são 1, com excessão de alpha_2
## que é igual a zero, pois é o expoente de B na interação
## ACD = A¹B°C¹D¹
## porque o efeito usado, ACD, não contém B. Portanto o vetor alpha fica
alpha <- c(1, 0, 1, 1)
db$bloco2 <- apply(db[, 1:4], 1, contr.def, alpha = alpha)
db
     A B C D bloco bloco2
(1)  0 0 0 0     0      0
a    1 0 0 0     1      1
b    0 1 0 0     1      0
ab   1 1 0 0     0      1
c    0 0 1 0     1      1
ac   1 0 1 0     0      0
bc   0 1 1 0     0      1
abc  1 1 1 0     1      0
d    0 0 0 1     1      1
ad   1 0 0 1     0      0
bd   0 1 0 1     0      1
abd  1 1 0 1     1      0
cd   0 0 1 1     0      0
acd  1 0 1 1     1      1
bcd  0 1 1 1     1      0
abcd 1 1 1 1     0      1
db[order(db$bloco2), ]
     A B C D bloco bloco2
(1)  0 0 0 0     0      0
b    0 1 0 0     1      0
ac   1 0 1 0     0      0
abc  1 1 1 0     1      0
ad   1 0 0 1     0      0
abd  1 1 0 1     1      0
cd   0 0 1 1     0      0
bcd  0 1 1 1     1      0
a    1 0 0 0     1      1
ab   1 1 0 0     0      1
c    0 0 1 0     1      1
bc   0 1 1 0     0      1
d    0 0 0 1     1      1
bd   0 1 0 1     0      1
acd  1 0 1 1     1      1
abcd 1 1 1 1     0      1

Exemplo

Exemplo 14-7, pg 357, Montgomery, Estatística aplicada e probabilidade para Engenheiros. Um experimento é realizado para investigar o efeito de quatro fatores sobre o desvio, em relação ao alvo, no disparo de um míssil. Os quatro fatores são: A = tipo de alvo, B = tipo de rastreador, C = altitude do alvo, D = distância do alvo. Cada fator pode ser convenientemente testado em 2 níveis e o sistema ótimo de rastreamento permitirá medir o desvio no disparo com a precisão de um pé. Dois atiradores diferentes são usados no teste de vôo e, já que há diferença entre operadores, os engenheiros de teste decidiram conduzir o planejamento \(2^4\) em 2 blocos com ABCD confundido. Faça a análise estatística dos dados.

##----------------------------------------------------------------------
## Resultados do experimento
da <- do.call(expand.grid, replicate(4, list(c(-1, 1))))
names(da) <- LETTERS[1:ncol(da)]
row.names(da) <- apply(da, 1,
                       function(i) paste(letters[1:4][i==1], collapse = ""))
row.names(da)[1] <- "(1)"
da$y <- c(3, 7, 5, 7, 6, 6, 8, 6, 4, 10, 4, 12, 8, 9, 7, 9)
da
      A  B  C  D  y
(1)  -1 -1 -1 -1  3
a     1 -1 -1 -1  7
b    -1  1 -1 -1  5
ab    1  1 -1 -1  7
c    -1 -1  1 -1  6
ac    1 -1  1 -1  6
bc   -1  1  1 -1  8
abc   1  1  1 -1  6
d    -1 -1 -1  1  4
ad    1 -1 -1  1 10
bd   -1  1 -1  1  4
abd   1  1 -1  1 12
cd   -1 -1  1  1  8
acd   1 -1  1  1  9
bcd  -1  1  1  1  7
abcd  1  1  1  1  9
##----------------------------------------------------------------------
## Definido os blocos

## Pela tabela de sinais
da$bloco <- with(da, A*B*C*D)

## Pelo contraste de definição
## Antes é necessário transformar a codificação para (0,1)
db <- as.data.frame(ifelse(da[, 1:4] == -1, 0, 1))
db$y <- da$y
db
     A B C D  y
(1)  0 0 0 0  3
a    1 0 0 0  7
b    0 1 0 0  5
ab   1 1 0 0  7
c    0 0 1 0  6
ac   1 0 1 0  6
bc   0 1 1 0  8
abc  1 1 1 0  6
d    0 0 0 1  4
ad   1 0 0 1 10
bd   0 1 0 1  4
abd  1 1 0 1 12
cd   0 0 1 1  8
acd  1 0 1 1  9
bcd  0 1 1 1  7
abcd 1 1 1 1  9
alpha <- c(1, 1, 1, 1)
db$bloco <- apply(db[, 1:4], 1, contr.def, alpha = alpha)

## Apenas para ilustração e verificação:
all.equal(row.names(da[order(da$bloco), ]),
          row.names(db[order(db$bloco, decreasing = TRUE), ]))
[1] TRUE
## Note que por qualquer um dos métodos a determinação de blocos fica a
## mesma. Daqui pra frente tanto faz usar um ou outro. O importante é
## identificar que a coluna de blocos deve ser um fator para poder
## entrar como um termo no modelo
da$bloco <- as.factor(da$bloco)

## Aqui procedemos da mesma forma. A diferença é que como estamos
## colocando bloco explicitamente no modelo, e bloco está confundido com
## a interação ABCD, então esta última interação não é especificada, e
## por isso, especificamos o modelo com todas as interações até terceira
## ordem apenas
X <- model.matrix(~ bloco + (A + B + C + D)^3, data = da)
colnames(X)
 [1] "(Intercept)" "bloco1"      "A"           "B"           "C"          
 [6] "D"           "A:B"         "A:C"         "A:D"         "B:C"        
[11] "B:D"         "C:D"         "A:B:C"       "A:B:D"       "A:C:D"      
[16] "B:C:D"      
## Calcula os contrastes, excluindo a interação e o bloco
contr <- t(X[, -(1:2)]) %*% da$y

## Efeitos = contraste/(n2^{k-1})
n <- 1 # sem repetições
k <- 4
ef <- contr/(n * 2^(k - 1))

## Gráfico de probabilidade normal dos efeitos
aux <- qqnorm(ef, col = 2, pch = 19); qqline(ef)
text(aux$x, aux$y, rownames(aux$y), cex = 0.8, pos = 3)

## Ajuste do modelo com interações de até segunda ordem
m0 <- lm(y ~ bloco + (A + B + C + D)^2, data = da)
anova(m0)
Analysis of Variance Table

Response: y
          Df  Sum Sq Mean Sq F value   Pr(>F)   
bloco      1  0.0625  0.0625  0.0588 0.820294   
A          1 27.5625 27.5625 25.9412 0.007016 **
B          1  1.5625  1.5625  1.4706 0.291974   
C          1  3.0625  3.0625  2.8824 0.164789   
D          1 14.0625 14.0625 13.2353 0.022003 * 
A:B        1  0.0625  0.0625  0.0588 0.820294   
A:C        1 22.5625 22.5625 21.2353 0.009969 **
A:D        1 10.5625 10.5625  9.9412 0.034416 * 
B:C        1  0.5625  0.5625  0.5294 0.507158   
B:D        1  0.5625  0.5625  0.5294 0.507158   
C:D        1  0.0625  0.0625  0.0588 0.820294   
Residuals  4  4.2500  1.0625                    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## A partir da ANOVA, vamos manter apenas os efeitos importantes
m1 <- update(m0, . ~ bloco + A + C + D + A:C + A:D)
anova(m1)
Analysis of Variance Table

Response: y
          Df  Sum Sq Mean Sq F value    Pr(>F)    
bloco      1  0.0625  0.0625  0.0796 0.7841601    
A          1 27.5625 27.5625 35.1239 0.0002217 ***
C          1  3.0625  3.0625  3.9027 0.0796325 .  
D          1 14.0625 14.0625 17.9204 0.0021961 ** 
A:C        1 22.5625 22.5625 28.7522 0.0004551 ***
A:D        1 10.5625 10.5625 13.4602 0.0051644 ** 
Residuals  9  7.0625  0.7847                      
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Análise dos resíduos
par(mfrow=c(2,2)); plot(m1); layout(1)

Confundimento do fatorial \(2^k\) em quatro blocos

É possível construir um fatorial \(2^k\) confundido em quatro blocos de \(2^{k-2}\) observações em cada. Estes experimentos são particularmente úteis quando o número de fatores é moderadamente alto (\(k \geq 4\)), e o tamanho dos blocos é relativamente pequeno.

Como exemplo considere um experimento \(2^5\). Se cada bloco deve conter 8 corridas, então são necessários 4 blocos. Para selecionar as combinações de tratamento em cada bloco, devemos selecionar dois efeitos a serem confundidos com blocos, por exemplo, ADE e BCE. Estes efeitos possuem os contrastes de definição

\[ \begin{align} L_1 = x_1 + x_4 + x_5 \\ L_2 = x_2 + x_3 + x_5 \end{align} \]

associados à eles. Com isso, cada combinação de tratamento irá gerar um particular par de valores de \(L_1 \mod 2\) e \(L_2 \mod 2\), ou seja, \((L_1, L_2) = (0,0), (0,1), (1,0), (1,1)\). Combinações de tratamentos que resultem no mesmo par de valores \((L_1, L_2)\) serão designadas para o mesmo bloco. No exemplo acima, portanto,

\[ \begin{align} L_1 = 0, L_2 = 0& \quad \text{para} \quad \text{(1), ad, bc, abcd, abe, ace, cde, bde} \\ L_1 = 1, L_2 = 0& \quad \text{para} \quad \text{a, d, abc, bcd, be, abde, ce, acde} \\ L_1 = 0, L_2 = 1& \quad \text{para} \quad \text{b, abd, c, acd, ae, de, abce, bcde} \\ L_1 = 1, L_2 = 1& \quad \text{para} \quad \text{e, ade, bce, abcde, ab, bd, ac, cd} \end{align} \]

Estas combinações de tratamentos serão atribuídas para cada um dos quatro blocos.

Devemos notar também que outra interação, além de ADE e BCE, deve estar confundida com bloco. Como existem 4 blocos e 3 graus de liberdade entre eles, e como ADE e BCE possuem 1 grau de liberdade cada, deve haver um efeito adicional com 1 grau de liberdade, que também deve estar confundido. Este efeito é a interação generalizada de ADE e BCE, que é definida como o produto de ADE com BCE. Portanto, no exemplo acima, a interação generalizada \((ADE)(BCE) = ABCDE^2 = ABCD\) (lembre-se das propriedades da tabela de sinais) também está confundida com blocos.

O procedimento geral para construir um experimento fatorial \(2^k\) em quatro blocos, é escolher dois efeitos para gerar os blocos, e automaticamente confundir um terceiro efeito que é a interação generalizada dos dois efeitos iniciais. Então o experimeno é construído utilizando-se os dois contrastes de definição \((L_1, L_2)\) para designar as combinações de tratamentos aos blocos.

Ao selecionar os efeitos a serem confudidos, devemos tomar cuidado para não escolher efeitos que sejam de interesse. Por exemplo, em um fatorial \(2^5\), podemos escolher confundir ABCDE e ABD, que automaticamente também confunde CE [\((ABCDE)(ABD) = A^2B^2CD^2E = CE\)], um efeito que provavelmente possa ser de interesse. Uma escolha melhor é confundir ADE com BCE, que automaticamente confunde ABCD. É preferível sacrificar informação de duas interações de terceira ordem, do que uma interação de segunda ordem.

##======================================================================
## Experimento fatorial 2^5 em quatro blocos
da <- do.call(expand.grid, replicate(5, list(c(0,1))))
names(da) <- LETTERS[1:ncol(da)]
row.names(da) <- apply(da, 1,
                       function(i) paste(letters[1:5][i==1], collapse = ""))
row.names(da)[1] <- "(1)"
da
      A B C D E
(1)   0 0 0 0 0
a     1 0 0 0 0
b     0 1 0 0 0
ab    1 1 0 0 0
c     0 0 1 0 0
ac    1 0 1 0 0
bc    0 1 1 0 0
abc   1 1 1 0 0
d     0 0 0 1 0
ad    1 0 0 1 0
bd    0 1 0 1 0
abd   1 1 0 1 0
cd    0 0 1 1 0
acd   1 0 1 1 0
bcd   0 1 1 1 0
abcd  1 1 1 1 0
e     0 0 0 0 1
ae    1 0 0 0 1
be    0 1 0 0 1
abe   1 1 0 0 1
ce    0 0 1 0 1
ace   1 0 1 0 1
bce   0 1 1 0 1
abce  1 1 1 0 1
de    0 0 0 1 1
ade   1 0 0 1 1
bde   0 1 0 1 1
abde  1 1 0 1 1
cde   0 0 1 1 1
acde  1 0 1 1 1
bcde  0 1 1 1 1
abcde 1 1 1 1 1
##----------------------------------------------------------------------
## Usando o contraste de definição, com ADE e BCE confundidos 

## L_1 = x1 + x_4 + x5
alpha1 <- c(1, 0, 0, 1, 1)
L1 <- apply(da, 1, contr.def, alpha = alpha1)
## L_2 = x2 + x_3 + x5
alpha2 <- c(0, 1, 1, 0, 1)
L2 <- apply(da, 1, contr.def, alpha = alpha2)

## Cria os blocos
da$bloco <- interaction(L1, L2, sep = "")
da <- da[order(da$bloco), ]
da
      A B C D E bloco
(1)   0 0 0 0 0    00
bc    0 1 1 0 0    00
ad    1 0 0 1 0    00
abcd  1 1 1 1 0    00
abe   1 1 0 0 1    00
ace   1 0 1 0 1    00
bde   0 1 0 1 1    00
cde   0 0 1 1 1    00
a     1 0 0 0 0    10
abc   1 1 1 0 0    10
d     0 0 0 1 0    10
bcd   0 1 1 1 0    10
be    0 1 0 0 1    10
ce    0 0 1 0 1    10
abde  1 1 0 1 1    10
acde  1 0 1 1 1    10
b     0 1 0 0 0    01
c     0 0 1 0 0    01
abd   1 1 0 1 0    01
acd   1 0 1 1 0    01
ae    1 0 0 0 1    01
abce  1 1 1 0 1    01
de    0 0 0 1 1    01
bcde  0 1 1 1 1    01
ab    1 1 0 0 0    11
ac    1 0 1 0 0    11
bd    0 1 0 1 0    11
cd    0 0 1 1 0    11
e     0 0 0 0 1    11
bce   0 1 1 0 1    11
ade   1 0 0 1 1    11
abcde 1 1 1 1 1    11
## Croqui do experimento
matrix(row.names(da), ncol = 4,
       dimnames = list(1:8, paste("Bloco", 1:4)))
  Bloco 1 Bloco 2 Bloco 3 Bloco 4
1 "(1)"   "a"     "b"     "ab"   
2 "bc"    "abc"   "c"     "ac"   
3 "ad"    "d"     "abd"   "bd"   
4 "abcd"  "bcd"   "acd"   "cd"   
5 "abe"   "be"    "ae"    "e"    
6 "ace"   "ce"    "abce"  "bce"  
7 "bde"   "abde"  "de"    "ade"  
8 "cde"   "acde"  "bcde"  "abcde"
##----------------------------------------------------------------------
## Também é possível fazer pela tabela de sinais
## Usando ADE e BCE cofundidos
db <- do.call(expand.grid, replicate(5, list(c(-1, 1))))
names(db) <- LETTERS[1:ncol(db)]
row.names(db) <- apply(db, 1,
                       function(i) paste(letters[1:5][i==1], collapse = ""))
row.names(db)[1] <- "(1)"
db
       A  B  C  D  E
(1)   -1 -1 -1 -1 -1
a      1 -1 -1 -1 -1
b     -1  1 -1 -1 -1
ab     1  1 -1 -1 -1
c     -1 -1  1 -1 -1
ac     1 -1  1 -1 -1
bc    -1  1  1 -1 -1
abc    1  1  1 -1 -1
d     -1 -1 -1  1 -1
ad     1 -1 -1  1 -1
bd    -1  1 -1  1 -1
abd    1  1 -1  1 -1
cd    -1 -1  1  1 -1
acd    1 -1  1  1 -1
bcd   -1  1  1  1 -1
abcd   1  1  1  1 -1
e     -1 -1 -1 -1  1
ae     1 -1 -1 -1  1
be    -1  1 -1 -1  1
abe    1  1 -1 -1  1
ce    -1 -1  1 -1  1
ace    1 -1  1 -1  1
bce   -1  1  1 -1  1
abce   1  1  1 -1  1
de    -1 -1 -1  1  1
ade    1 -1 -1  1  1
bde   -1  1 -1  1  1
abde   1  1 -1  1  1
cde   -1 -1  1  1  1
acde   1 -1  1  1  1
bcde  -1  1  1  1  1
abcde  1  1  1  1  1
## Obtém as intreações de confundimento
db$ADE <- with(db, A*D*E)
db$BCE <- with(db, B*C*E)

## Cria os blocos
db$bloco <- with(db, interaction(ADE, BCE, sep = ""))
db <- db[order(db$bloco), ]
db
       A  B  C  D  E ADE BCE bloco
(1)   -1 -1 -1 -1 -1  -1  -1  -1-1
bc    -1  1  1 -1 -1  -1  -1  -1-1
ad     1 -1 -1  1 -1  -1  -1  -1-1
abcd   1  1  1  1 -1  -1  -1  -1-1
abe    1  1 -1 -1  1  -1  -1  -1-1
ace    1 -1  1 -1  1  -1  -1  -1-1
bde   -1  1 -1  1  1  -1  -1  -1-1
cde   -1 -1  1  1  1  -1  -1  -1-1
a      1 -1 -1 -1 -1   1  -1   1-1
abc    1  1  1 -1 -1   1  -1   1-1
d     -1 -1 -1  1 -1   1  -1   1-1
bcd   -1  1  1  1 -1   1  -1   1-1
be    -1  1 -1 -1  1   1  -1   1-1
ce    -1 -1  1 -1  1   1  -1   1-1
abde   1  1 -1  1  1   1  -1   1-1
acde   1 -1  1  1  1   1  -1   1-1
b     -1  1 -1 -1 -1  -1   1   -11
c     -1 -1  1 -1 -1  -1   1   -11
abd    1  1 -1  1 -1  -1   1   -11
acd    1 -1  1  1 -1  -1   1   -11
ae     1 -1 -1 -1  1  -1   1   -11
abce   1  1  1 -1  1  -1   1   -11
de    -1 -1 -1  1  1  -1   1   -11
bcde  -1  1  1  1  1  -1   1   -11
ab     1  1 -1 -1 -1   1   1    11
ac     1 -1  1 -1 -1   1   1    11
bd    -1  1 -1  1 -1   1   1    11
cd    -1 -1  1  1 -1   1   1    11
e     -1 -1 -1 -1  1   1   1    11
bce   -1  1  1 -1  1   1   1    11
ade    1 -1 -1  1  1   1   1    11
abcde  1  1  1  1  1   1   1    11
## Veja que o resultado é o mesmo
cbind(row.names(da), row.names(db))
      [,1]    [,2]   
 [1,] "(1)"   "(1)"  
 [2,] "bc"    "bc"   
 [3,] "ad"    "ad"   
 [4,] "abcd"  "abcd" 
 [5,] "abe"   "abe"  
 [6,] "ace"   "ace"  
 [7,] "bde"   "bde"  
 [8,] "cde"   "cde"  
 [9,] "a"     "a"    
[10,] "abc"   "abc"  
[11,] "d"     "d"    
[12,] "bcd"   "bcd"  
[13,] "be"    "be"   
[14,] "ce"    "ce"   
[15,] "abde"  "abde" 
[16,] "acde"  "acde" 
[17,] "b"     "b"    
[18,] "c"     "c"    
[19,] "abd"   "abd"  
[20,] "acd"   "acd"  
[21,] "ae"    "ae"   
[22,] "abce"  "abce" 
[23,] "de"    "de"   
[24,] "bcde"  "bcde" 
[25,] "ab"    "ab"   
[26,] "ac"    "ac"   
[27,] "bd"    "bd"   
[28,] "cd"    "cd"   
[29,] "e"     "e"    
[30,] "bce"   "bce"  
[31,] "ade"   "ade"  
[32,] "abcde" "abcde"
##======================================================================
## Experimento fatorial 2^4 em quatro blocos
da <- do.call(expand.grid, replicate(4, list(c(0,1))))
names(da) <- LETTERS[1:ncol(da)]
row.names(da) <- apply(da, 1,
                       function(i) paste(letters[1:4][i==1], collapse = ""))
row.names(da)[1] <- "(1)"
da
     A B C D
(1)  0 0 0 0
a    1 0 0 0
b    0 1 0 0
ab   1 1 0 0
c    0 0 1 0
ac   1 0 1 0
bc   0 1 1 0
abc  1 1 1 0
d    0 0 0 1
ad   1 0 0 1
bd   0 1 0 1
abd  1 1 0 1
cd   0 0 1 1
acd  1 0 1 1
bcd  0 1 1 1
abcd 1 1 1 1
##----------------------------------------------------------------------
## Usando o contraste de definição, com ABC e ACD confundidos, por
## consequência, a interação generalizada é
## (ABC)(ACD) = A²BC²D = BD

## L_1 = x1 + x_2 + x3
alpha1 <- c(1, 1, 1, 0)
L1 <- apply(da, 1, contr.def, alpha = alpha1)
## L_2 = x1 + x_3 + x4
alpha2 <- c(1, 0, 1, 1)
L2 <- apply(da, 1, contr.def, alpha = alpha2)

## Cria os blocos
da$bloco <- interaction(L1, L2, sep = "")
da <- da[order(da$bloco), ]
da
     A B C D bloco
(1)  0 0 0 0    00
ac   1 0 1 0    00
abd  1 1 0 1    00
bcd  0 1 1 1    00
b    0 1 0 0    10
abc  1 1 1 0    10
ad   1 0 0 1    10
cd   0 0 1 1    10
ab   1 1 0 0    01
bc   0 1 1 0    01
d    0 0 0 1    01
acd  1 0 1 1    01
a    1 0 0 0    11
c    0 0 1 0    11
bd   0 1 0 1    11
abcd 1 1 1 1    11
## Croqui do experimento
matrix(row.names(da), ncol = 4,
       dimnames = list(1:4, paste("Bloco", 1:4)))
  Bloco 1 Bloco 2 Bloco 3 Bloco 4
1 "(1)"   "b"     "ab"    "a"    
2 "ac"    "abc"   "bc"    "c"    
3 "abd"   "ad"    "d"     "bd"   
4 "bcd"   "cd"    "acd"   "abcd" 

Exemplos

##======================================================================
## Montgomery, Design and analysis of experiments, 5 ed., ex. 7-7, com
## dados do ex. 6-21
url <- "http://www.leg.ufpr.br/~fernandomayer/dados/montgomery_6-21.txt"
dados <- read.table(url, header = TRUE)
row.names(dados) <- apply(dados, 1,
                       function(i) paste(letters[1:5][i==1], collapse = ""))
row.names(dados)[1] <- "(1)"
dados
       A  B  C  D  E  y
(1)   -1 -1 -1 -1 -1  7
a      1 -1 -1 -1 -1  9
b     -1  1 -1 -1 -1 34
ab     1  1 -1 -1 -1 55
c     -1 -1  1 -1 -1 16
ac     1 -1  1 -1 -1 20
bc    -1  1  1 -1 -1 40
abc    1  1  1 -1 -1 60
d     -1 -1 -1  1 -1  8
ad     1 -1 -1  1 -1 10
bd    -1  1 -1  1 -1 32
abd    1  1 -1  1 -1 50
cd    -1 -1  1  1 -1 18
acd    1 -1  1  1 -1 21
bcd   -1  1  1  1 -1 44
abcd   1  1  1  1 -1 61
e     -1 -1 -1 -1  1  8
ae     1 -1 -1 -1  1 12
be    -1  1 -1 -1  1 35
abe    1  1 -1 -1  1 52
ce    -1 -1  1 -1  1 15
ace    1 -1  1 -1  1 22
bce   -1  1  1 -1  1 45
abce   1  1  1 -1  1 65
de    -1 -1 -1  1  1  6
ade    1 -1 -1  1  1 10
bde   -1  1 -1  1  1 30
abde   1  1 -1  1  1 53
cde   -1 -1  1  1  1 15
acde   1 -1  1  1  1 20
bcde  -1  1  1  1  1 41
abcde  1  1  1  1  1 63
##======================================================================
## Analisando os dados originais, SEM considerar blocos ainda

## Grafico de probabilidade normal dos efeitos
tab <- model.matrix(~ A*B*C*D*E, data = dados)
contr <- t(tab[, -1]) %*% dados$y
n <- 1
k <- 5
ef <- contr/(n * 2^(k-1))
aux <- qqnorm(ef, col = 2, pch = 19); qqline(ef)
text(aux$x, aux$y, rownames(aux$y), cex = 0.8, pos = 3)

## ANOVA
## Efeitos importantes: A, B, C, e AB
m0 <- lm(y ~ A + B + C + A:B, data = dados)
anova(m0)
Analysis of Variance Table

Response: y
          Df Sum Sq Mean Sq F value    Pr(>F)    
A          1 1116.3  1116.3  382.27 < 2.2e-16 ***
B          1 9214.0  9214.0 3155.34 < 2.2e-16 ***
C          1  750.8   750.8  257.10 2.534e-15 ***
A:B        1  504.0   504.0  172.61 3.038e-13 ***
Residuals 27   78.8     2.9                      
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Analise dos residuos
qqnorm(residuals(m0)); qqline(residuals(m0))

##======================================================================
## Ex. 7-7. Dividindo o experimento em 2 blocos, usando ABCD como
## confundimento

## Usando o contraste de definição
## Antes, converte para a codificação (0,1)
dados[, 1:5] <- ifelse(dados[, 1:5] == -1, 0, 1)
dados
      A B C D E  y
(1)   0 0 0 0 0  7
a     1 0 0 0 0  9
b     0 1 0 0 0 34
ab    1 1 0 0 0 55
c     0 0 1 0 0 16
ac    1 0 1 0 0 20
bc    0 1 1 0 0 40
abc   1 1 1 0 0 60
d     0 0 0 1 0  8
ad    1 0 0 1 0 10
bd    0 1 0 1 0 32
abd   1 1 0 1 0 50
cd    0 0 1 1 0 18
acd   1 0 1 1 0 21
bcd   0 1 1 1 0 44
abcd  1 1 1 1 0 61
e     0 0 0 0 1  8
ae    1 0 0 0 1 12
be    0 1 0 0 1 35
abe   1 1 0 0 1 52
ce    0 0 1 0 1 15
ace   1 0 1 0 1 22
bce   0 1 1 0 1 45
abce  1 1 1 0 1 65
de    0 0 0 1 1  6
ade   1 0 0 1 1 10
bde   0 1 0 1 1 30
abde  1 1 0 1 1 53
cde   0 0 1 1 1 15
acde  1 0 1 1 1 20
bcde  0 1 1 1 1 41
abcde 1 1 1 1 1 63
## Como o fator a ser confundido é ABCD, então
## L = x_1 + x_2 + x_3 + x_4 + x_5
alpha <- c(1, 1, 1, 1, 1)
dados$bloco <- apply(dados[, 1:5], 1, contr.def, alpha = alpha)
dados <- dados[order(dados$bloco), ]
dados
      A B C D E  y bloco
(1)   0 0 0 0 0  7     0
ab    1 1 0 0 0 55     0
ac    1 0 1 0 0 20     0
bc    0 1 1 0 0 40     0
ad    1 0 0 1 0 10     0
bd    0 1 0 1 0 32     0
cd    0 0 1 1 0 18     0
abcd  1 1 1 1 0 61     0
ae    1 0 0 0 1 12     0
be    0 1 0 0 1 35     0
ce    0 0 1 0 1 15     0
abce  1 1 1 0 1 65     0
de    0 0 0 1 1  6     0
abde  1 1 0 1 1 53     0
acde  1 0 1 1 1 20     0
bcde  0 1 1 1 1 41     0
a     1 0 0 0 0  9     1
b     0 1 0 0 0 34     1
c     0 0 1 0 0 16     1
abc   1 1 1 0 0 60     1
d     0 0 0 1 0  8     1
abd   1 1 0 1 0 50     1
acd   1 0 1 1 0 21     1
bcd   0 1 1 1 0 44     1
e     0 0 0 0 1  8     1
abe   1 1 0 0 1 52     1
ace   1 0 1 0 1 22     1
bce   0 1 1 0 1 45     1
ade   1 0 0 1 1 10     1
bde   0 1 0 1 1 30     1
cde   0 0 1 1 1 15     1
abcde 1 1 1 1 1 63     1
## Croqui do experimento
matrix(row.names(dados), ncol = 2,
       dimnames = list(1:16, paste("Bloco", 1:2)))
   Bloco 1 Bloco 2
1  "(1)"   "a"    
2  "ab"    "b"    
3  "ac"    "c"    
4  "bc"    "abc"  
5  "ad"    "d"    
6  "bd"    "abd"  
7  "cd"    "acd"  
8  "abcd"  "bcd"  
9  "ae"    "e"    
10 "be"    "abe"  
11 "ce"    "ace"  
12 "abce"  "bce"  
13 "de"    "ade"  
14 "abde"  "bde"  
15 "acde"  "cde"  
16 "bcde"  "abcde"
## Analise dos dados
## Note que para calcular os contrastes, precisamos dos dados na
## codificação (-1,1)
dados[, 1:5] <- ifelse(dados[, 1:5] == 0, -1, 1)
## Monta a tebela de sinais. Repare que as interações vão até a quata
## ordem, pois a interação de quinta ordem, ABCDE, está confundida com
## bloco, que já está presente no modelo
tab <- model.matrix(~ bloco + (A+B+C+D+E)^4, data = dados)
contr <- t(tab[, -(1:2)]) %*% dados$y
n <- 1
k <- 5
ef <- contr/(n * 2^(k-1))
aux <- qqnorm(ef, col = 2, pch = 19); qqline(ef)
text(aux$x, aux$y, rownames(aux$y), cex = 0.8, pos = 3)

## Efeitos importantes: A, B, C, e AB
m1 <- lm(y ~ bloco + A + B + C + A:B, data = dados)
anova(m1)
Analysis of Variance Table

Response: y
          Df Sum Sq Mean Sq   F value    Pr(>F)    
bloco      1    0.3     0.3    0.0931    0.7627    
A          1 1116.3  1116.3  369.4296 < 2.2e-16 ***
B          1 9214.0  9214.0 3049.3532 < 2.2e-16 ***
C          1  750.8   750.8  248.4686 8.027e-15 ***
A:B        1  504.0   504.0  166.8075 8.080e-13 ***
Residuals 26   78.6     3.0                        
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Analise dos residuos
qqnorm(residuals(m1)); qqline(residuals(m1))