5.1. Análise de Variância Multivariada em Blocos Completos Casualizados
Em experimentos cujo delineamento utilizado é o de blocos completos casualizados (EBCC), e cada parcela é avaliada sob as respostas de duas ou mais variáveis, os vetores de médias dos tratamentos devem ser comparados, usando-se a técnica da Análise de Variância Multivariada (MANOVA).
Considere-se o quadro genérico de observações de um experimento multivariado planejado em blocos completos casualizados.
Quadro 1. Valores genéricos dos dados observados em um EBCC
Tratamento 1 Tratamento t Totais de Blocos
Bloco V1 ... Vp ... V1 ... Vp V1 ... Vp Bloco 1 y111 ... y1p1 ... yt11 ... ytp1 B11 ... Bp1 Bloco2 y112 ... y1p2 ... yt12 ... ytp2 B12 ... Bp2 .... ... ... ... ... ... ... ... ... ... ... Bloco b y11b ... y1pb ... yt1b ... ytpb B1b ... Bpb Totais de Tratamento T11 ... T1p ... Tt1 ... Ttp Médias de Tratamentos y11. ... y1p. ... yt1. ... ytp. 5.2. Modelo 5.2.1. Modelo Linear
A representação dos dados observacionais (Quadro 1), obtidos a partir de experimentos agronômicos distribuídos em blocos completos casualizados com p respostas exógenas pode ser descrita pelo modelo probabilístico com estrutura de erros dependentes para uma dada resposta, como:
y
ijk= µ
j+ β
jk+ τ
ij+ ε
ijk ondej = 1, 2, ..., p é o índice da variável resposta; k = 1, 2, ..., b é o índice de bloco;
y
ijk é o valor observado da j-ésima variável resposta, do i-ésimo tratamento no k- ésimo bloco;µ
j é o efeito comum da variável resposta Vj;β
jk é o efeito do k-ésimo bloco na j-ésima variável resposta;τ
ijé o efeito do i-ésimo tratamento na variável resposta Vj eε
ijk é o erro aleatório, relativo ao valor observadoy
ijk. O vetor de errosε
ik tem distribuição normal multivariada, com vetor nulo de médias e matriz comum de covariância Σ (p x p), quaisquer que sejam o tratamento e o bloco considerados. Portanto,ε'
ik = (ε
i1k, ε
i2k, ..., ε
ipk) ~ Np (0, Σ)para qualquer i = 1, 2, ..., t e k = 1, 2, ..., b. Além disso os vetores
ε'
ik, correspondentes às diferentes unidades amostrais em qualquer bloco, são independentemente distribuídos.5.2.2 Modelo matricial
O modelo matricial dos dados descritos anteriormente pode ser indicado como:
no qual
Y é a matriz das observações, de dimensões (tb x p), onde as b primeiras linhas representam os elementos do tratamento 1, assim sendo até as b últimas linhas correspondentes ao tratamento t; y111 y121 ... y1p1 y112 y122 ... y1p2 y11b y12b ... y1pb y211 y221 ... y2p1 y212 y222 ... y2p2 Y = y21b y22b ... y2pb yt11 yt21 ... ytp1 yt12 yt22 ... ytp2 yt1b yt2b .... ytpb (tb x p);
A é a matriz de delineamento do experimento, formada por uns e zeros, tipo [tb x (b+t+1), sendo as b primeiras colunas correspondentes ao parâmetro β, tendo 1 no bloco em consideração e zero no restante; a coluna seguinte corresponde ao parâmetro τ, tendo 1 no tratamento especificado e zero nos demais e, finalmente, a última coluna constituída de uns, correspondendo ao parâmetro µ.
: : : : : : : :
1 0 ... 0 1 0 ... 0 1 0 1 ... 0 1 0 ... 0 1 A = . 1 0 ... 0 0 0 ... 1 1 0 1 … 0 0 0 ... 1 1 0 0 .... 1 0 0 .... 1 1 tb x (b+t+1)
ξ é a matriz de parâmetros populacionais, do tipo (b+t+1) x p, e pode ser descrita da seguinte forma: : : : : : :
β11 β12 .... β1p
βb1 βb2 .... βbp
ξ = τ11 τ12 .... τ1p
τt1 τt2 .... τtp
µ1 µ2 .... µp (b + t + 1) x p
ε
é a matriz, de ordem bt x p, que contém as variações residuais e111 e121 ... e1p1 e112 e122 ... e1p2 e11b e12b ... e1pb e211 e221 ... e2p1 e212 e222 ... e2p2ε
= e21b e22b ... e2pb et11 et21 ... etp1 et12 et22 ... etp2 et1b et2b .... etpb (tb x p) : : : : : : : : : : : :5.3. MANOVA
Na análise de variância multivariada em blocos completos casualizados foi realizado um teste de hipóteses geral, referente ao efeito de tratamento. A hipótese de nulidade do teste pode ser descrita como:
H0: não existe efeito de tratamento
ou, equivalentemente:
τ11 τ21 τt1
: : :
H0: : = : = ... = :
τ1p τ2p τtp
Matricialmente, a hipótese testada pode ser descrita da seguinte forma: H0: C ξ U = 0
no qual,
C de ordem [(t-1) x (b+t+1)], é a matriz que realiza operações fundamentais entre os níveis de tratamento e pode ser escrita de diferentes maneiras, entre elas:
0 .... 0 1 -1 0 .... 0 0 0
0 .... 0 0 1 -1 .... 0 0 0
C (t-1) x (b+t+1)
=
...
e U é a matriz, que à semelhança de C, envolve as variáveis respostas nas suas operações fundamentais. No presente trabalho, a matriz U será a unidade de ordem p.
Considere uma partição apropriada da matriz C= [C1 C2] no qual as
submatrizes têm as dimensões[(t-1)x(b+t)] e [(t-1)x1], respectivamente. Da mesma forma, uma partição semelhante na matriz de parâmetros ξ, ou seja, ξ = [ξ1 ξ2]. Sob estas
condições, a hipótese equivalente a ser testada pode ser reescrita por: H0: C1 ξ1 U + C2 ξ2 U = 0
Partição equivalente deve ser feita em A, estabelecendo novos conjuntos de matrizes com partições: [C1 C2], [ξ1 ξ2] e [A1 A2], cujas formas genéricas
encontram-se no Apêndice 1.
O teste de hipóteses do efeito de tratamentos pode ser realizado considerando as raízes características de E-1H, onde as matrizes H e E são dados por:
H = U’Y’A1(A1’A1)-1C1[C1(A1’A1)-1C1]-1(A1’A1)-1A1’YU
E = U’Y’[I - A1(A1’A1)-1 A1]YU
As matrizes H e E correspondem, respectivamente, às somas de quadrados e produtos de Tratamento e Resíduo, e seus elementos podem ser representados aleatoriamente pelas expressões usuais do cálculo de soma de quadrados e produtos dos valores observados.
Através da utilização da Álgebra de Matrizes, têm-se os elementos genéricos das matrizes simétricas H = [hrs], B = [brs], T = [trs] e E = [ers], todas de
ordem (pxp) dados por:
hrs = (1/b) Σ Tir Tis – 1/(bt) Gr Gs
;
i=1brs = (1/t) Σ Brk Bsk – 1/(bt) Gr Gs
;
trs = Σ Σ yirk yisk – 1/(bt) Gr Gs
;
ers= trs – hrs – brs, onde r,s = 1, 2, ..., p.
no qual,
Bjk = Σ yijk (total do bloco k relativo à variável Vj
;
j = 1, ..., p e k = 1, ..., b);Tij = Σ yijk (total do tratamento i para Vj; j = 1, ..., p e i = 1, ..., t);
Gj = Σ Σ yijk
=
Σ Bjk=
Σ Tij (total geral relativo à variável Vj, j = 1, ..., p).Uma forma mais comum de expressar as matrizes H, B, T e E é designar os elementos genéricos através de somas de quadrados e produtos, obtidos de forma idêntica às descritas anteriormente.
Ou seja:
SQTot (1) SPTot (1,2) SPTot (1,3) ... SPTot (1,p) SQTot (2) SPTot (2,3) ... SPTot (2,p) T(p x p) = SQTot (3) ... SPTot (3,p) ...…... SQTot (p) i=1 k=1 t t s k=1 b k=1 b i=1 t k=1 b k=1 b i=1 t i=1 t
SQTrat (1) SPTrat (1,2) SPTrat (1,3) ... SPTrat (1,p) SQTrat (2) SPTrat (2,3) ... SPTrat (2,p) H(p x p) = SQTrat (3) ... SPTrat (3,p) ...……… SQTrat (p) SQBl (1) SPBl (1,2) SPBl (1,3) ... SPBl (1,p) SQBl (2) SPBl (2,3) ... SPBl (2,p) B (p x p) = SQBl (3) ... SPBl (3,p) ...…….. SQBl (p)
SQRes (1) SPRes (1,2) SPRes (1,3) ... SPRes (1,p) SQRes (2) SPRes (2,3) ... SPRes (2,p)
E
(p x p)=
SQResl (3) ... SPRes (3,p)... SQRes (p)
As matrizes envolvidas no teste estatístico referente ao efeito de tratamento podem ser dispostas no Quadro da Análise de Variância Multivariada (MANOVA), como segue.
Quadro 2. Quadro da Análise de Variância Multivariada em EBCC
Causa de variação Matriz
Tratamentos H
Blocos B
Resíduo E
Total T
no qual:
T é a matriz de Soma de Quadrados e Produtos Total
H é a matriz de Soma de Quadrados e Produtos de Tratamento B é a matriz de Soma de Quadrados e Produtos de Bloco E é a matriz de Soma de Quadrados e Produtos Resíduo.
Para a construção da estatística do teste de hipóteses, considera-se a equação polinomial |H - λE| = 0, ou equivalentemente, |H E-1 - λI| = 0
As raízes da equação polinomial são chamadas de raízes características de H E-1, as quais no presente estudo serão indicadas por: λ
1 ≥ λ2 ≥ ... λs > 0 (raízes
características não nulas).
A partir das raízes características, várias estatísticas, sob a veracidade de H0, podem ser definidas para o teste de hipóteses:
Com os parâmetros:
m = (|p – q| -1)/2 ; n = [t (b-1)-b-p] / 2 e s = min (p, q) = número de raízes características não nulas. q = t – 1;
ii) Critério de Wilks: Λ = Π (1 - θl);
iii) Traço de Lawley-Hotelling: T2 = Σ λl ;
iv) Traço de Pillai: V = Σθl ;
onde λl são as s raízes características não-nulas de H E-1. Deve ser
destacado que o teste de Roy é um procedimento exato.
GODOI (1985) realiza as seguintes aproximações para um valor correspondente da estatística de F, a qual possibilita usar a tabela de Snedecor/Fisher para a regra de decisão.
Para o Critério de Wilks é considerado, neste trabalho, W = [(m1 s1 – 2a) (1 - Λ1/s1)] / [p (t – 1) Λ1/s1)]; sendo que
para p ≥ 3 e q ≥ 3, W é aproximadamente distribuído como F[(p – 1) (b – 1) , ||m1s1 – 2a ||]
onde m1 = b (t – 1) – [½ (p + t)];
a = [(p – 1) (b – 1) - 2]/ 4;
s1 =
{
[p2 (b – 1)2 - 4] / [p2 +(t – 1)2 - 5]}
1/2 e||m1 s1 – 2 a||, onde || . || indica o maior inteiro que não supera ao seu valor.
MORRISON (1981) considera, para os experimentos que apresentam exatamente duas variáveis respostas (bivariada), a seguinte estatística do critério de Wilks, utilizado na MANOVA: l=1 s l=1 s l=1 s
Λ = det E / det ( H + E )
Para avaliar a significância do resultado de Λ, é realizada a seguinte aproximação:
F = [(1 - Λ1/2) / (Λ1/2)] [(2n + 2) / (2m + 3)]
com distribuição F de Fischer/Snedecor, com (4m + 6) e 4(n + 1) graus de liberdade, sob a veracidade de H0.
Da mesma forma, tem-se para o Critério do Traço de Hotelling- Lawley:
HL = [2 (s n + 1) T2] / [s2 (2 m + s + 1)], cuja aproximação é realizada pela distribuição F s (2m + s + 1), 2 (sn + 1).
Na mesma linha de procedimento, o Traço de Pillai considera: P = [(2 n + s + 1) V] / [(2 m + s + 1) (s – V)]
sendo que P é aproximadamente distribuído como F s (2m + s + 1), s (2n + s + 1).
Para o uso das aproximações de Morrison os parâmetros m, n e s; assim como as estatísticas Λ, T2 e V são as definidas anteriormente.
Em todas as considerações, a hipótese de ausência de efeito de tratamento pode ser testada, bastando para isso comparar a estatística aproximada com o valor crítico (nível de significância) tabelado Fα e utilizar a regra usual de decisão, ou seja, não se
rejeitar H0 se F aproximado ≤ F tabelado e, rejeitar H0, se F calculado > F tabelado. Quando se
concretiza a rejeição da hipótese da nulidade, torna-se interessante a construção de intervalos de confiança para contrastes entre os efeitos de tratamentos, com nível de confiança de 100 (1 - α)%.
Para a construção dos limites (inferior e superior) do intervalo com 100 (1 - α)% de confiança, para contrastes entre pares de médias de dois tratamentos, em uma dada variável resposta, considere:
i) a′j = (0 ... 1 ... 0) um vetor p dimensional, de elementos nulos, exceção a j-ésima
posição que será preenchida com a unidade.
ii) Os vetores de médias dos tratamentos yi, para os contrastes são:
yi..’ = (yi1. y12. ... y1p.)
para
yij.= (1/b) Σ yijk
i = 1, 2, ..., t .
iii) Os limites de confiança de Roy. O primeiro método adotado neste trabalho para a construção do intervalo de confiança para contrastes entre pares de tratamentos é o princípio de união-intersecção de Roy, sendo que os limites, inferior e superior, são calculados, respectivamente:
LI: aj′ (yi. - yi’.) - {[ (Xα / (1 – Xα) ] aj′ E aj (2/b) }1/2;
LS: aj′ (yi. - yi..) + {[ (Xα / (1 – Xα) ] aj′ E aj (2/b) }1/2;
j = 1, ..., p e i , i’ = 1, …, t; no qual
Xα é o valor crítico de Roy, no nível α de significância, podendo ser obtido dos ábacos de Heck
ou das tabelas de Pillai, com parâmetros (Apêndice 4):
m = (|p – q| -1)/2 ; n = [t (b-1)-b-p] / 2 ; s = min (p, q) = número de raízes características não nulas e
k=1 b
q = t – 1.
Os limites de confiança referem-se ao contraste entre o i-ésimo e i′-ésimo tratamentos, fixada a j-ésima variável resposta.
iv) O segundo método adotado neste trabalho para o estudo de intervalo de confiança é a estatística de Bonferroni e os limites são obtidos através das seguintes fórmulas:
LI: aj′ (yi. - yi’.) – t {ϕres ; α/[tp (t -1)]} [aj′ E aj (2/b) (1/ϕres)] 1/2
LS: aj′ (yi. - yi’.) + t{ϕres ; α/[tp (t -1)]} [aj′ E aj (2/b) (1/ϕres)] 1/2
no qual:
ϕres = (b-1)(t-1)
t{ϕres ; α/[tp (t -1)]} é o valor correspondente ao quantil de ordem 100 {1 – [α/tp(t-1)]}% da
distribuição t de Student com ϕres graus de liberdade (Apêndice 5).
Variando-se o vetor aj de maneira apropriada, em ambos os casos,
poderão ser construídos todos os intervalos de interesse para a diferença entre duas médias. Para a interpretação da significância do contraste deve-se proceder da seguinte maneira: se o intervalo de confiança excluir o zero, o contraste estudado será considerado significativo, no nível α pré-estabelecido; caso contrário, não há razão estatística que identifique diferença entre as médias dos tratamentos, na variável em discussão.
5.4. Programa Computacional – MANOVASYS - EBCC
Na elaboração do programa computacional para a metodologia utilizou-se o Delphi, uma ferramenta para programação visual estruturada baseada em Object Pascal, que por sua vez é uma linguagem muito robusta e fortemente indicada para desenvolvimento de sistemas de grande porte. A plataforma utilizada no desenvolvimento foi a Windows, tornando assim, o sistema compatível com as plataformas Windows 95, Windows 98 e Windows ME, que são os sistemas operacionais mais utilizados no momento, e de aplicação em diversas áreas científicas.
O programa apresenta opção de entrada de dados via-teclado ou arquivo-texto armazenado em disquete ou gravado no disco rígido. Os dados para importação podem estar digitados na planilha Excel ou em outro editor e salvos em texto formatado, separado por espaços (prn). No primeiro caso, os dados deverão ser digitados fazendo uso da vírgula; no segundo caso, o programa realizará a conversão do ponto em vírgula, caso o usuário tenha utilizado o ponto na planilha Excel.
Uma tela com os valores digitados se abrirá e, ao se clicar o botão “Prosseguir”, uma Lista com os resultados parciais surgirá fornecendo os parâmetros necessários para identificar os valores críticos, necessários para o cálculo dos Intervalos de Confiança do Princípio da União-Intersecção de Roy e de Bonferroni e com os resultados das estatísticas. Se pelo menos uma das estatísticas for significativa, será realizado o cálculo de F para cada uma das variáveis isoladamente. O usuário é dispensado do fornecimento dos valores críticos da distribuição probabilística F utilizada, quando trabalhar no nível comum de significância (5%), pois esses valores encontram-se embutido no programa, sendo que, se
optar por outros níveis de interesse, o usuário necessitará utilizar a tabela de valores da distribuição probabilística.
Para a obtenção dos Intervalos de Confiança do Princípio da União- Intersecção de Roy, é necessário fornecer os valores críticos, Xα e t, respectivamente. Esses
valores podem ser obtidos no Apêndice ou no próprio programa.
Ao final do cálculo dos intervalos de confiança, o usuário poderá salvar toda análise na Pasta do próprio ManovaSys-EBCC.
O procedimento para a utilização do programa MANOVASYS - EBCC encontra-se detalhado no Manual do Usuário no Apêndice 2.
5.5. Esquema Experimental do Exemplo Ilustrativo
Os dados utilizados para aplicação do método são adaptados de DEMÉTRIO (1995) e referem-se a experimentos realizados com cana-de-açúcar por técnicos do CTC (Centro Tecnológico Copersucar), Piracicaba-SP, tendo as seguintes características:
Delineamento : Blocos Completos Casualizados Número de Variedades (Tratamentos): 4
Variedades utilizadas : 1. IAC 48-65 2. IAC 51-205 3. IAC 52-150 4. IAC 53-37 Número de áreas (Blocos) : 8
Unidade experimental : 5 linhas de 8 metros
Os ensaios foram conduzidos em dois locais diferentes, de acordo com o Quadro 3:
Quadro 3. Locais de realização dos experimentos de cana-de-açúcar Local
Usina Fazenda
Instalação Corte Colheita
São José Z. L. Patos 07/02/80 1º 01/07/81
As variáveis analisadas, segundo interesse da pesquisadora, são: fibra% (V1), brix% (V2), pol% (V3) e produtividade (t/ha) (V4). Para a obtenção dos dados das características tecnológicas, é tomada uma amostra de 10 canas por parcela, a qual é passada por prensa hidráulica, dando origem ao bolo úmido e ao caldo. O bolo úmido, pesado e levado à estufa, dá origem ao bolo seco, obtendo-se, por razão, a porcentagem de fibra da cana, ou seja,
fibra% de cana = [(Peso do bolo seco) / (Peso do bolo úmido)] * 100%.
A partir do caldo, por medição direta com refratômetro digital, é obtida a variável brix%, e após a clarificação do caldo é obtida a variável pol % do caldo, com sacarímetro automático e, então,
pol % de cana = P (1 – 0,01 F) 0,945 % onde,
P = pol % do caldo; F = fibra % de cana.
A variável produtividade (produção em t/ha) é obtida colhendo-se e pesando-se a parcela toda.