#========================================================================================== # Controle de Processos Industriais www.leg.ufpr.br/ce074 # Curso de Estatística - 2012/2 # Prof. Walmes Zeviani # LEG - DEST - UFPR # # Aula 14 (08/02/2013) # Experimentos fatoriais 3^k, técnicas de blocagem e fracionamento. #========================================================================================== #------------------------------------------------------------------------------------------ # Representação dos efeitos com "1 GL". ## Cada fator agora tem 3 níveis e com isso temos 2 graus de liberdade em cada fator. ## Representa-se esses 2 graus de liberdade em contrastes que são ortogonais entre sí ## e que para fatores de níveis quantitativos e equidistantes representam efeito linear ## e quadrático. O contraste linear tem [-1, 0, 1] e quadrático [1, -2, 1] como ## coeficientes. Efeitos e somas de quadrados são calculadas da mesma forma que no 2^k. ## Como cada fator tem 2 GL, as interações duplas têm 2*2=4 GL, as triplas 2*2*2=8 e assim ## por diante. As interações duplas podem ser decompostas em 4 contrastes ortogonais vindos ## dos produtos dos efeitos L e Q dos fatores envolvidos. ## da <- expand.grid(A1=-1:1, B1=-1:1) # efeitos linear ## da <- transform(da, A2=3*A1^2-2, B2=3*B1^2-2) # efeitos quadráticos ## X <- model.matrix(~(A1+A2)*(B1+B2), data=da) # matriz com o produto dos termos ## X # A com 2 GL, B com 2 GL e AB com 4 GL #------------------------------------------------------------------------------------------ # Cálculo das somas de quadrados. da$y <- rnorm(nrow(da)) # um vetor de valores observados qualquer ## sq = contraste^2/soma_pesos^2 contr <- t(X)%*%da$y den <- apply(X, 2, function(x) sum(x^2)) # equivalente à diag(crossprod(X)) sq <- contr^2/den #------------------------------------------------------------------------------------------ # Representação "2 GL", útil para técnicas de confundimento. ## Os níveis dos fatores são representados por 0, 1, e 2. As "palavras" agora terão ## expoentes 0, 1 e 2 nas letras que representam os fatores. No fatorial 2^k, devido ter ## apenas e níveis esses expoentes eram 0 e 1 e nem se percebia a sua presença. As ## "palavras" terão a seguinte forma A^pB^q para o fatorial 3^2. Para a primeira letra ## sempre o expoente será 1. Assim para o fatorial 3^2 representamos os 4 de liberdade ## da interação AB por dois conjuntos de contrastes ortogonais A¹B¹ e A¹B², ou omitindo ## o expoente 1, AB e AB², ambos com 2 GL cada. Em um 3³ temos AB = A¹B¹ e A¹B², ## AC = A¹C¹ e A¹C² e BC = B¹C¹ e B¹C². ABC = A¹B¹C¹, A¹B²C¹, A¹B¹C² e A¹B²C², cada um ## com 2 GL. Esse tipo de representação parece complicada mas é extremamente útil para ## fazer blocagem e frações. #------------------------------------------------------------------------------------------ # Blocagem no experimento 3^k, em 3 blocos. ## Assim como no 2^k temos que especificar um efeito para ser confundido com o efeito de ## blocos. Como serão usados 3 blocos, temos 2 GL envolvidos neles. Os efeitos com a ## representação acima (de expoentes) também têm 2 GL. Aí que está lógica, vamos confundir ## um desses efeitos com o de blocos, tal como é no 2^k. Considere um 3² = 9 corridas. ## Temos os seguintes efeitos que podem ser estimados: A, B, A¹B¹ e A¹B², cada um deles ## com 2 GL, perfazendo 9-1=8. O mais adequado é confunfir uma interação alta, no caso ## A¹B² será considerada para formar o contraste de definição ## L = 1*x1+2*x2+0*x3 (mod 3). ## Como mod 3 só retorna os valores 0, 1 e 2, esses são os falores usados para atribuir ## as corridas aos blocos. da <- expand.grid(A=0:2, B=0:2) rl <- function(x, l) sum(x*l)%%3 da$bloco <- apply(da, 1, rl, l=c(1,2)) da[order(da$bloco),] db <- expand.grid(A=0:2, B=0:2, C=0:2) # relação de definição A¹B²C² db$bloco <- apply(db, 1, rl, l=c(1,2,2)) split(db, db$bloco) #------------------------------------------------------------------------------------------ # Fatorial fracionário no experimento 3^k ## Segue a mesma regra de seleção de corridas usada na blocagem. ## Para determinar os efeitos confundidos é semelhante ao 2^k. Pegamos a relação de ## definição e fazermos o "produto" pelo efeito que desejamor determinar o confundimento. ## Por exemplo, no 3³ com I = A¹B²C², para verificar quem está confundido com A fazemos ## A*I = A*(A¹B²C²) = A²B²C² = (A²B²C²)² = A^4 B^4 C^4 = A¹B¹C¹, lembre que 4 mod 3 = 1; ## A*I = A*(A¹B²C²)² = A³B^4 C^4 = BC = B¹C¹ ## A interação AB ## AB*I = AB*(A¹B²C²) = A²B³C² = (A²B³C²)² = A^4 B^6 C^4 = A¹C¹ ## AB*I = AB*(A¹B²C²)² = A³B^5 C^4 = (A³B^5 C^4)² = A^6 B^10 C^8 = B¹C² ## A regra é elevar a expressão ao quadrado o quanto for necessário até que o expoente de ## A fique 0 ou 1. #------------------------------------------------------------------------------------------ # Análise de um 3³ dividido em 3 blocos, I = A¹B²C². da <- read.table("http://www.leg.ufpr.br/~walmes/data/pimentel_npk2.txt", header=TRUE, sep="\t") str(da) da$bloco <- with(da, interaction(bloco, rept)) str(da) require(lattice) xyplot(producao~N|K, groups=P, data=da, type=c("p","a")) xtabs(~K+N+P, data=da) m0 <- lm(producao~bloco+poly(N,2)*poly(P,2)*poly(K,2), data=da) anova(m0) #------------------------------------------------------------------------------------------