Técnicas de diagnóstico de influência local na análise espacial da produtividade da soja

Diagnostic techniques of local influence in spatial analysis of soybean yield

ARTIGOS CIENTÍFICOS

SANEAMENTO E CONTROLE AMBIENTAL

Técnicas de diagnóstico de influência local na análise espacial da produtividade da soja

Diagnostic techniques of local influence in spatial analysis of soybean yield

Joelmir A. Borssoi^I; Miguel A. Uribe-Opazo^II; Manuel Galea^III

^IInstituto de Matemática e Estatística, Universidade de São Paulo, São Paulo- SP, [e-mail]

^IICentro de Ciências Exatas e Tecnológicas, Universidade Estadual do Oeste do Paraná, Cascavel - PR

^IIIDepartamento de Estatística, Pontifícia Universidad Católica de Chile, Santiago, Chile

RESUMO

A modelagem da estrutura de dependência espacial pela abordagem da geoestatística é fundamental para a definição de parâmetros que definem esta estrutura, e que são utilizados na interpolação de valores em locais não amostrados pela técnica de krigagem. Entretanto, a estimação de parâmetros pode ser muito afetada pela presença de observações atípicas nos dados amostrados. O desenvolvimento deste trabalho teve por objetivo utilizar técnicas de diagnóstico de influência local em modelos espaciais lineares gaussianos, utilizados em geoestatística, para avaliar a sensibilidade dos estimadores de máxima verossimilhança e máxima verossimilhança restrita na presença de dados discrepantes. Estudos com dados experimentais mostraram que tanto a presença de valores atípicos como de valores considerados influentes, pela análise de diagnóstico, pode exercer forte influência nos mapas temáticos, alterando, assim, a estrutura de dependência espacial. As aplicações de técnicas de diagnóstico de influência local devem fazer parte de toda análise geoestatística a fim de garantir que as informações contidas nos mapas temáticos tenham maior qualidade e possam ser utilizadas com maior segurança pelo agricultor.

Palavras-chave: geoestatística, máxima verossimilhança, máxima verossimilhança restrita

ABSTRACT

Modeling of spatial dependence structure, concerning geoestatistics approach, is an indispensable tool for fixing parameters that define this structure, applied on interpolation of values in places that are not sampled, by kriging techniques. However, the estimation of parameters can be greatly affected by the presence of atypical observations on sampled data. Thus, this trial aimed at using diagnostics techniques of local influence in spatial linear Gaussians models, applied at geoestatistics in order to evaluate sensitivity of maximum likelihood estimators and restrict maximum likelihood to small perturbations in these data. So, studies with simulated and experimental data were performed. Those results, obtained from the study of real data, allowed us to conclude that the presence of atypical values among the sampled data can have a strong influence on thematic maps, changing, therefore, the spatial dependence. The application of diagnostics techniques of local influence should be part of any geoestatistic analysis, ensuring that the information from thematic maps has better quality and can be used with greater security by farmers.

Keywords: geoestatistic, maximum likelihood, restrict maximum likelihood

INTRODUÇÃO

Nas últimas décadas, vem aumentando a busca por tecnologias que melhorem a qualidade das informações obtidas com respeito a uma área agrícola seja pela busca de melhores equipamentos, seja pela busca de novas formas de análise de dados.

A geoestatística, que tem como fundamento a teoria das variáveis regionalizadas, é um método que considera a distribuição espacial das medidas, permitindo definir o raio de auto-correlação espacial entre elementos amostrais e, consequentemente, a distância máxima onde se tem dependência espacial. O estudo da variabilidade espacial é considerado de grande importância, segundo CRESSIE (1993), e não considerar a variabilidade espacial pode impedir que diferenças reais sejam levantadas, o que impossibilitaria um tratamento diferenciado de acordo com as exigências locais.

Na modelagem da variabilidade espacial, são estimados os parâmetros que definem a estrutura de dependência espacial e que serão utilizados na interpolação de valores, em locais não amostrados, pela técnica de krigagem. Assim, são gerados mapas temáticos que poderão ser utilizados para uma aplicação localizada de insumos ou manejo localizado do solo. A qualidade destes mapas depende da qualidade das inferências dos modelos ajustados. Mas, para que a interpolação produza previsões confiáveis e represente a real variabilidade local, o processo de modelagem deve ser feito com muita cautela, principalmente na presença de pontos discrepantes ou influentes.

Nesta perspectiva, o objetivo deste trabalho foi utilizar técnicas de diagnóstico de influência local em modelos espaciais lineares gaussianos para verificar a existência de observações que possam exercer algum tipo de influência nas inferências sobre os parâmetros que definem a estrutura de dependência espacial e na construção dos mapas temáticos. Para tanto, realizaram-se estudos de influência local, utilizando os métodos de estimação de parâmetros de máxima verossimilhança (MV) e de máxima verossimilhança restrita (MVR) a um conjunto de dados da produtividade da soja, tendo como covariáveis a resistência mecânica do solo à penetração e a densidade do solo em diferentes profundidades.

MATERIAL E MÉTODOS

Modelos espaciais lineares

Para modelar dados com estrutura espacial, considere-se, segundo MARDIA & MARSHALL (1984), um processo estocástico gaussiano {Z(s), s∈S}, com , sendo o espaço euclidiano, d-dimensional (d > 1). Suponha-se que os dados, Z(s₁),...,Z(s_n), deste processo sejam registrados em localizações espaciais conhecidas, s_i (i = 1,...,n) , e gerados pelo seguinte modelo:

em que, os termos determinístico µ(s_i) e estocástico ε(s_i) podem depender da localização espacial onde Z(s_i) foi obtida.

Assume-se que o erro estocástico ε (.) tem média zero, E[ε (s_i)] = 0, e que a variação entre pontos no espaço é determinada por alguma função de covariância C(s_i,s_u) = COV[ε (s_i), ε (s_u)] e que, para algumas funções conhecidas de s, como x₁(s),..., x_p(s), a média do processo estocástico é:

em que,

β₁,...,β_p - parâmetros desconhecidos e a serem estimados.

De forma equivalente, em notação matricial, tem-se:

Então,

E(ε)=0(vetor nulo) e a matriz de covariância de Σ = [(σ_iu)], onde σ_iu = C(s_i , s_u). Assume-se que Σ é não singular, que X tem colunas com posto completo e que Z segue uma distribuição gaussiana multivariada com média Xβ e matriz de covariância Σ, isto é, Z ~ N_n (Xβ, Σ).

Considerando-se de maneira particular a forma paramétrica da matriz de covariância, segundo

em que,

φ₁ - efeito pepita, ou erro de variância;

φ₂ - contribuição, ou variância de dispersão;

R - matriz que é função de φ₃, R = R(φ₃) = [(r_iu)], matriz n x n simétrica com seus elementos da diagonal r_ii = 1, i = 1, ...,n, em que φ₃ é função do alcance (a) do modelo, e

I_n - matriz identidade n x n.

A forma paramétrica da matriz de covariância, dada na Equação 4, ocorre para vários processos isotrópicos, onde a covariância C(s_i, s_u) é definida segundo a função de covariâncias C(h_iu) = φ₂r_iu, em que h_iu = ||s_i – s_u|| é a distância euclidiana entre os pontos s_ie s_u. Nas funções de covariâncias C(h_iu), a variância do processo estocástico Z é C(0) = φ₁ + φ₂, e a semivariância pode ser definida como γ(h) = C(0) – C(h).

Influência local

Em muitas situações, pode-se observar um conjunto de dados com observações aberrantes ou discrepantes que podem ser consideradas influentes, isto é, que podem mudar algum tipo de decisão na construção de modelos geoestatísticos.

O método de influência local proposto por COOK (1986) avalia o efeito simultâneo de observações sobre os estimadores de MV e MVR sem a necessidade de sua eliminação do conjunto de dados. As técnicas de influência local têm-se constituído num procedimento útil para realizar análises de sensibilidade de modelos estatísticos e têm sido amplamente utilizadas em modelos de regressão linear e não linear. CHRISTENSEN et al. (1993) estudaram métodos de diagnóstico baseado em eliminação de casos em modelos lineares espaciais. BORSSOI et al. (2009) apresentam o gráfico de influência local |L_max| versus i (ordem dos dados) como medida de diagnóstico de influência local para a detecção de pontos influentes.

Para um conjunto de dados observados, seja l(θ) o logaritmo da função verossimilhança do modelo postulado, em que θ =(β^T, φ^T)^T, e seja ω um vetor de perturbações pertencente a um espaço de perturbações Ω. Seja l(θ/ω) o logaritmo da função verossimilhança correspondente ao modelo perturbado por ω ∈ Ω. Assume-se que existe um vetor ω₀∈ Ω tal que l(θ) = l(θ/ω₀), para todo θ e que l(θ/ω) é duas vezes diferenciável em (θ^T, ω^T)^T.

O logaritmo da função verossimilhança perturbada l(θ/ω), considerando a perturbação Z_ω = Z + ω para o modelo normal, é dado por

A influência da perturbação ω, nos estimadores MV e MVR do vetor de parâmetros θ, pode ser avaliada pelo afastamento da verossimilhança, definido por:

em que,

- estimador de máxima verossimilhança de θ do modelo postulado, e

- estimador de máxima verossimilhança de θ do modelo perturbado.

COOK (1986) propôs estudar o comportamento local de LD(ω) em torno de ω₀, utilizando a curvatura normal C_l de LD(ω) em ω₀, na direção de algum vetor unitário l, definido como sendo:

em que,

L - matriz de informação observada, avaliada em , e

Δ - matriz (p + q) x n dada por Δ = (Δ_β^T, Δ_φ^T)^T, avaliada em e em ω = ω_0, em que:

A matriz de informação L é definida como: , em que,

, com seus elementos:

Considere-se a matriz B = Δ^TL^-1Δ, e a i-ésima observação definida como

em que, b_ii são os elementos da diagonal principal da matriz B.

Pode-se utilizar o gráfico de C_i versus i (ordem dos dados) como técnica de influência local para avaliar a existência de observações influentes. Sugere-se considerar a i-ésima observação como influente se seu valor for maior que duas vezes a média das medidas C_i.

Para avaliar a influência local no preditor linear p(s₀, θ), considere-se que (s₀, θ) é um vetor n x 1:

Na eq.(10), tem-se que ,

em que, , com seus elementos dados por:

para j = 1,..., q.

A direção de máxima inclinação local é obtida pela eq.(11):

Pode-se utilizar a medida ou simplesmente l_p para analisar a sensibilidade do preditor linear.

Estudo experimental

Foram analisados os dados da variável Produtividade da soja (Prod) [t ha^-1] e as covariáveis: Densidade do solo (Des) [mg m^-3] e Resistência mecânica do solo à penetração (RSP) [MPa], nas camadas de 0 - 0,10 m, 0,10 - 0,20 m e 0,20 - 0,30 m de profundidade, nos pontos onde foi coletada a produtividade, obtidos no ano agrícola de 2004/2005, em uma área comercial de produção de grãos de 57 ha. A técnica de amostragem utilizada foi a sistemática alinhada, com distância de 75 m entre pontos.

A área está localizada no município de Cascavel, na região oeste do Paraná, aproximadamente 24,95º sul de latitude, 53,57º oeste de longitude e altitude média de 650 m. O solo é classificado como Latossolo Vermelho distroférrico, e o clima da região é temperado mesotérmico e superúmido, tipo climático Cfa (Köeppen), com temperatura anual média de 21 °C. Todos os elementos amostrais foram georreferenciados com auxílio de um receptor de sinal do GPS (Global Positioning System) marca Trimble, modelo GeoExplorer3, num sistema espacial de coordenadas (UTM).

Neste trabalho, utilizaram-se os modelos exponencial, gaussiano e família Matérn, com os métodos de estimação de parâmetros máxima verossimilhança (MV) e máxima verossimilhança restrita (MVR) (MARDIA & MASHALL, 1984). Para escolher o modelo espacial que melhor se ajusta à semivariância, utilizou-se a técnica de validação cruzada (Faraco et al., 2008) e o máximo valor do logaritmo da função verossimilhança (MLL). Para investigar a existência de pontos influentes, realizou-se uma análise de diagnóstico por meio da técnica de influência local. Para a análise dos dados, foi utilizado o software R (R DEVELOPMENT CORE TEAM, 2009) e o módulo geoestatístico geoR (RIBEIRO JR & DIGGLE, 2009).

RESULTADOS E DISCUSSÃO

Na Tabela 1, são apresentadas as estatísticas descritivas para a variável Produtividade da soja (Prod), Resistência mecânica do solo à penetração nas camadas de 0,0 a 0,10 m, 0,1 a 0,20 m e 0,20 a 0,3 m de profundidade (RSP_0,0-0,1, RSP_0,1-0,2 e RSP_0,2-0,3, respectivamente) e Densidade do solo nas camadas de 0,0 a 0,10 m, 0,10 a 0,20m e 0,20 a 0,30 m de profundidade (Des_0-0,1, Des_0,1-0,2 e Des_0,2-0,3, respectivamente). Verifica-se que o desvio-padrão (DP = 0,38 t ha^-1) da produtividade da soja é pequeno comparado com o valor da média, que foi de 3,23 t ha^-1, apresentando homogeneidade nos dados da produtividade pelo baixo valor do coeficiente de variação (C.V. < 30%). Observa-se também que existe homogeneidade nos dados de Des_0,0-0,1, Des_0,1-0,2, Des_0,2-0,3, RSP_0,0-0,1, RSP_0,1-0,2 e RSP_0,2-0,3 também pelo baixo valor do coeficiente de variação. Pelo gráfico de pontos espaciais classificados por quartis (Box-plot) da produtividade da soja (Figura 1a), observa-se a presença de um ponto discrepante, com valor igual a 2,09 t ha^-1, sendo o 13º valor da série de dados localizado na parte inferior esquerda (Figura 1b).

Análises de diagnóstico

Foram aplicadas técnicas gráficas de diagnóstico de influência local com a finalidade de avaliar, na produtividade da soja, se o ponto discrepante (13º valor) ou alguns outros estariam exercendo algum tipo de influência na estrutura de dependência espacial (afastamento da verossimilhança) e no preditor linear. As técnicas utilizadas foram apresentadas por meio dos gráficos dos coeficientes C_ie l_p obtidos das Equações (8) e (11).

Nas Figuras 2 e 3, são apresentadas as análises de influência local para a variável produtividade, com o total de 47 dados coletados (Prod), utilizando os modelos espaciais ajustados: exponencial, gaussiano e Matérn (com parâmetro kappa igual a 0,7) com os métodos de estimação de parâmetros MV e MVR. Na Tabela 2, é apresentado um resumo dos resultados das técnicas de diagnóstico de influência local com as respectivas observações consideradas influentes pelas Figuras 2 e 3, para todas as combinações entre os modelos exponencial, gaussiano e Matérn, e os métodos de estimação de parâmetros MV e MVR.

Pela Tabela 2, verifica-se que a observação número 13, apontada como discrepante pelo gráfico Box-Plot (Figura, 1a), foi identificada como influente no afastamento da verossimilhança em todos os casos (técnica C_i), exceto para modelo Matérn (k = 0,7) associado ao método de estimação MVR. Já na avaliação da influência no preditor linear (técnica l_p), ela foi identificada apenas pelo modelo exponencial associado ao método de estimação MV. As observações 3º, 8º, 13º, 23º e 30º foram identificadas com maior frequência na Tabela 2, por isso foram retiradas do conjunto total de dados e novamente realizadas as análises para avaliar suas influências. As observações retiradas da amostra de dados e seus valores são: 3º (3,64 t ha^-1), 8º (3,00 t ha^-1), 13º (2,09 t ha^-1), 23º (3,02 t ha^-1) e 30º (2,87 t ha^-1). Por simplicidade, definem-se as variáveis: Prod (amostra da produtividade e covariáveis com todas as observações coletadas) e Prod^* (amostra da produtividade e covariáveis sem as observações 3º, 8º, 13º, 23º e 30º).

a) Influência nas medidas descritivas

Pela Tabela 3, observa-se que há pouca variação nas medidas descritivas das variáveis Prod e Prod^*. As maiores diferenças são verificadas no valor mínimo (Mín.) e no coeficiente de variação (CV), em que seu valor teve uma redução, passando de 11,79% (Prod) para 10,62% (Prod^*), ambos indicando homogeneidade dos dados amostrais.

b) Influência nas estimativas dos parâmetros da modelagem espacial da produtividade da soja como função de covariáveis

JOHANN et al. (2004) concluíram que, dentre os atributos físicos dos solos estudados, a resistência mecânica do solo à penetração, na camada de 0,0 a 0,1 m de profundidade, foi a variável que melhor se correlacionou com a produtividade da soja. Assim, da definição de Z(s_i) na eq.(1), e do termo determinístico μ (s_i) na eq.(2), a média do processo estocástico é representada por uma relação linear com as covariáveis RSP e Des por meio da eq.(12):

em que,

β₁,..., β₇ - parâmetros desconhecidos e a serem estimados.

Nas Tabelas 4 e 5, são apresentados os resultados das análises de variabilidade espacial da produtividade da soja, por meio da geoestatística, para os dados originais (Prod) e retirando-se os pontos considerados influentes (Prod*).

Pelos resultados da Tabela 4, pode-se observar, para todos os modelos ajustados, quando retiradas as observações 3º, 8º, 13º, 23º e 30º (Prod^*), influenciaram na estimação dos parâmetros de β₂, β₄, β₅, β₆ e β₇, pois, quando estas foram retiradas (Prod^*), os valores das estimativas dos parâmetros aumentaram e o valor de β₃ diminuiu, com relação aos obtidos para Prod.

Observando a Tabela 5, verifica-se que, quando o parâmetro φ₁ (efeito pepita) foi estimado sem as observações influentes (Prod^*), em todos os modelos ajustados, houve redução no seu valor, quando comparado com os obtidos para a variável Prod. Um comportamento exatamente ao contrário verifica-se para o parâmetro φ₃ (função do alcance), mostrando que, na estimação de parâmetros, houve uma relação inversa entre φ₁ e φ₃ quando as observações influentes foram retiradas. Quanto à φ₂ (contribuição), dentro de todos os modelos, identifica-se o seguinte comportamento: quando a estimativa é feita por MV, o valor de φ₂ para Prod^* é menor do que o obtido para Prod; já por MVR, o valor de φ₂ para Prod^* é sempre maior do que o obtido para Prod.

As estimativas para o parâmetro φ₃ para Prod e Prod^* apresentam as maiores diferenças, tanto pelo método de MV quanto pelo método de MVR.

c) Influência na validação dos modelos

Na Tabela 6, são apresentados os resultados da validação de modelos ajustados para as variáveis Prod e Prod^*, utilizando o critério da validação cruzada (FARACO et al., 2008) e o valor máximo do logaritmo da função verossimilhança (MLL).

Para a variável Prod, foi escolhido como melhor ajuste o modelo Matérn (k=0,7) com parâmetros estimados por MVR; e quando foram retiradas as observações consideradas influentes pelas técnicas de influência local, para a variável Prod^*, foi escolhido como melhor ajuste o modelo gaussiano, com parâmetros estimados por MVR, utilizando os critérios da validação cruzada, e tendo o indicador MLL entre os valores máximos (Tabela 6). Isto mostra que as observações 3º, 8º, 13º, 23º e 30º são influentes também na escolha do melhor modelo espacial ajustado.

d) Influência na construção dos mapas temáticos

Na Figura 4, apresentam-se os mapas temáticos das variáveis Prod (Figura 4a) e Prod^* (Figura 4b), confeccionados por meio da interpolação por krigagem universal, visto que, neste estudo, foram introduzidas covariáveis. Os mapas foram construídos utilizando-se dos modelos indicados pelo critério da validação cruzada, que é o modelo Matérn com k= 0,7 para Prod e o modelo gaussiano para Prod* (Tabela 6).

Pela Figura 4, observa-se que há grande alteração na distribuição da produtividade da área em estudo, principalmente nas localizações das sub-regiões com produtividade nos intervalos de 3,04 a 3,23 t ha^-1 e de 3,23 a 3,43 t ha^-1.

Estas constatações indicam que as observações 3º, 8º, 13º, 23º e 30º estão exercendo grande influência também na construção dos mapas temáticos, confirmando o que a medida de diagnóstico no preditor linear (gráficos l_p) já havia indicado (Figura 3).

CONCLUSÕES

O estudo mostrou que a presença de valores atípicos entre dados amostrados pode exercer forte influência nos mapas temáticos, alterando, assim, a dependência espacial. Assim, o mapa da distribuição espacial da produtividade que seria apresentado ao agricultor, não estaria representando de forma correta as tendências da produtividade na área em estudo. Isso acarretaria em interpretações errôneas no gerenciamento da área quanto às futuras intervenções no trato do solo de forma localizada.

Conclui-se, também, que nem sempre um valor discrepante é influente, ao passo que um valor qualquer, mesmo não discrepante, pode ser considerado influente. Isto é, um valor discrepante não deve ser eliminado da amostra antes de ser realizada uma análise de diagnóstico de influência local.

A aplicação de técnicas de diagnóstico de influência local deve fazer parte de toda análise geoestatística, garantindo que as informações contidas nos mapas temáticos tenham maior qualidade e possam ser utilizadas com maior segurança pelo agricultor.

AGRADECIMENTOS

Ao CNPq, CONAB e Fundação Araucária, pelo apoio financeiro, e à UNIOESTE - Universidade Estadual do Oeste do Paraná.

Recebido pelo Conselho Editorial em: 12-9-2009

Aprovado pelo Conselho Editorial em: 25-2-2011

• CHRISTENSEN, R.; JOHNSON, W.; PEARSON, L. Covariance function diagnostics for spatial linear models. International Association for Mathematical Geology, New York, v.25, n.2, p.145-160, 1993.
• BORSSOI, J.A.; URIBE-OPAZO, M.A.; GALEA, M. Diagnostic techniques applied in geostatistics for agricultural data analysis. Revista Brasileira de Cięncia do Solo, Viçosa-MG, v.33, n.6, p. 1.561-1.570, 2009.
• COOK, R.D. Assessment of local influence (with discussion). Journal of the Royal Statistical Society, Series B, London, v.48, n.2,  p.133-169, 1986.
• CRESSIE, N.A.C. Statistic for spatial data New York: John Wiley & Sons, 1993. 900 p.
• FARACO M.A.; URIBE-OPAZO, M.A.; SILVA, E.A.; JOHANN J.A.; BORSSOI, J.A. Seleçăo de modelos de variabilidade espacial para elaboraçăo de mapas temáticos de atributos físicos do solo e produtividade da soja. Revista Brasileira de Cięncia do Solo, Viçosa-MG, v.32, n.2, p.463-476, 2008.
• JOHANN, J.A.; URIBE-OPAZO, M.A.; SOUZA, E.G.; ROCHA, J.V. Variabilidade espacial dos atributos físicos do solo e da produtividade em uma área experimental de agricultura de precisăo com e sem manejo localizado. Revista Brasileira de Engenharia Agrícola e Ambiental, Campina Grande, v.8, n.2, p.1-8. 2004.
• MARDIA, K.; MARSHALL, R. Maximum likelihood estimation of models for residual covariance in spatial regression. Biometrika, Alemanha, v.71, p.135-146, 1984.
• R DEVELOPMENT CORE TEAM R: A language and environment for statistical computing. R Foundation for statistical Computing, Vienna, Austria. ISBN 3-90005107-0. Disponível em: <http://www.r-project.org>. Acesso em: 20 jul. 2009
• RIBEIRO JR., P.J.; DIGGLE P. J. geoR: A package for geostatistical analysis. R-NEWS, 01. Disponível em: <http://cran.r-project.org/doc/Rnews 2001>. Acesso em: 3 jun. 2009.

Datas de Publicação

• Publicação nesta coleção
  06 Maio 2011
• Data do Fascículo
  Abr 2011

Histórico

• Recebido
  12 Set 2009
• Aceito
  25 Fev 2011