%«echo=F»= %require(geoR) %rm(list=ls(all=TRUE)) %@

\subsection{Análise da variável alumínio medida na camada~1}

% % \subsubsection{Análise exploratória dos dados} % %

\noindent Criando o objeto geodata para Alumínio. Estamos indicando que as coordenadas $(x,y)$ dos pontos amostrados encontram-se na primeira e segunda coluna do banco de dados (arquivo Quimicos.txt), a variável Alumínio encontra-se na terceira coluna e a covariável região encontra-se na coluna de número~12.

«»= Alum←as.geodata(Quimicos,coords.col=1:2,data.col=3,covar.col=12) @

\noindent O comando a seguir atribui à variável {\sf Alum\$borders} as coordenadas que se encontram no arquivo {\sf borda.txt}.

«»= Alum$borders ← borda @ %$

\noindent Descrevendo os dados rapidamente «»= summary(Alum) @

\noindent Visualizando graficamente os dados originais, ou seja, sem qualquer transformação ou arranjo. «fig=T»= plot(Alum) @

\noindent Pelo histograma podemos observar que os dados não seguem uma distribuição normal, sendo necessária uma transformação nos dados. Para efetuarmos essa transformação, precisaremos identificar o valor de $\lambda$ e aplicar o método de Box-Cox. A figura \ref{fig:boxcoxAl} dá uma indicação gráfica de quais valores de $\lambda$ pertencem ao intervalo de 95\% de confiança. No comando seguinte, por conveniência visual, o gráfico foi construído em um intervalo de -1 a 0.5. O resultado foi armazenado na variável {\sf BC1} (sem o efeito da região) e {\sf BC2} (com o efeito da região).

«fig=F»= par(mfrow=c(1,2)) boxcox(Alum,lam = seq(-1,0.5, l=100))→ BC1 boxcox(Alum, trend=~regiao, lam = seq(-1,0.5, l=100))→ BC2 par(mfrow=c(1,1)) @

\begin{figure}[!ht] \SweaveOpts{width=7,height=4} \centering «fig=T, echo=F»= par(mfrow=c(1,2)) boxcox(Alum,lam = seq(-1,0.5, l=100))→ BC1 boxcox(Alum, trend=~regiao, lam = seq(-1,0.5, l=100))→ BC2 par(mfrow=c(1,1)) @ \vspace{-0.5cm} \caption{Log-verossimilhança do parâmetro $\lambda$ para a transformação de Box-Cox (esquerda) sem efeito da região e (direita) com efeito da região} \label{fig:boxcoxAl} \end{figure}

\noindent Considerando um nível de 99\% de confiança decidimos adotar $\lambda=0$ devido a proximidade dos valores originais a zero, usando a log-normal e evitando assim problemas no truncamento.

\noindent Revendo, graficamente, os dados incluindo um efeito da região temos:

\begin{figure}[!ht] \SweaveOpts{width=5.5,height=5.5} «fig=T»= plot(Alum, trend=~regiao) @ \end{figure}

% % % \subsubsection{Análise Espacial} % % \noindent Inspecionando os dados espacialmente (Figura \ref{fig:postplotAl} ). Neste ponto, consideramos o valor de $\lambda=0$ e estudaremos o efeito da região no modelo.

«fig=F»= par(mar=c(3.5,3,0,0.5),mgp=c(2,0.8,0)) points(Alum, lambda=0,cex.min=.3,cex.max=1, axes=F, ann=F,col=gray(seq(0.9,0,l=length(Alum$data)))) polygon(A1) polygon(A2) par(mar=c(3.5,3.5,0.5,0.5),mgp=c(2,0.8,0)) @ %$

«fig=F»= par(mar=c(3.5,3,0,0.5),mgp=c(2,0.8,0)) points(Alum, trend=~regiao, lambda=0,cex.min=.3,cex.max=1, axes=F, ann=F,col=gray(seq(0.9,0,l=length(Alum$data)))) polygon(A1) polygon(A2) par(mar=c(3.5,3.5,0.5,0.5),mgp=c(2,0.8,0)) @ %$

\begin{figure}[!ht] \centering \SweaveOpts{width=6,height=3.5} \setkeys{Gin}{width=0.90\textwidth} «fig=T,echo=FALSE»= par(mfrow=c(1,2), mar=c(1,2,1,1),mgp=c(2,0.8,0)) points(Alum, lambda=0,cex.min=.3,cex.max=1, axes=F, ann=F,col=gray(seq(0.9,0,l=length(Alum$data)))) polygon(A1) polygon(A2)

points(Alum, trend=~regiao, lambda=0,cex.min=.3,cex.max=1, axes=F, ann=F,col=gray(seq(0.9,0,l=length(Alum$data)))) polygon(A1) polygon(A2) par(mfrow=c(1,2), mar=c(2,2,1,1),mgp=c(2,0.8,0)) @ %$ \vspace{-1.5cm} \caption{Locação dos dados na área de estudo e suas respectivas regiões, sem e com a covariável} \label{fig:postplotAl} \end{figure}

\noindent A próxima ação é a construção do semivariograma experimental a partir dos dados transformados, visando uma análise exploratória do comportamento da correlação espacial. Foi considerado também o efeito das diferentes regiões. A figura \ref{fig:variogramaAl} sugere uma dependência espacial, com valores iniciais de efeito pepita de 0,6, variância total de 1,4 e um alcance de 150 metros.

«fig=F, results=hide»= Alum.v←variog(Alum,max.dist=350,lambda=0, trend=~regiao) plot(Alum.v, pch=19,cex=0.8) @ \begin{figure}[!ht] \centering \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.60\textwidth} «fig=T,echo=FALSE, results=hide»= Alum.v←variog(Alum,max.dist=350,lambda=0, trend=~regiao) plot(Alum.v, pch=19,cex=0.8) @ \vspace{-.5cm} \caption{Semivariograma experimental com agrupamento de semivariâncias em intervalos de distâncias} \label{fig:variogramaAl} \end{figure}

% % \subsubsection{Ajustando um modelo teórico ao semivariograma experimental} % %

\noindent A idéia aqui é avaliar a aderência de uma função de correlação de Matérn ao semivariograma experimental produzido pelos dados transformados, segundo 5 valores diferentes do parâmetro de suavidade {\it kappa}. Os parâmetros aqui são estimados pelo método da máxima-verossimilhança:

\begin{itemize} \item sem o efeito da covariável {\sf região}. «results=hide»= Al.ml ← list() Al.ml$Al.1←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=0.5, lambda=0) Al.ml$Al.2←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=1.5, lambda=0) Al.ml$Al.3←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=2.5, lambda=0) Al.ml$Al.4←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=3.5, lambda=0) Al.ml$Al.5←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=4.5, lambda=0) @ %$

\item com o efeito da covariável {\sf região}. «results=hide»= Al.ml1 ← list() Al.ml1$Al.1←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=0.5, lambda=0, trend=~regiao) Al.ml1$Al.2←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=1.5, lambda=0, trend=~regiao) Al.ml1$Al.3←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=2.5, lambda=0, trend=~regiao) Al.ml1$Al.4←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=3.5, lambda=0, trend=~regiao) Al.ml1$Al.5←likfit(Alum, ini=c(0.6,150),cov.model="matern", kappa=4.5, lambda=0, trend=~regiao) @ %$ \end{itemize}

\noindent Assim, podemos apresentar os resultados obtidos da otimização da função log-verossimilhanca e escolher o modelo que apresentar o maior valor para a log-verossimilhança.

\begin{itemize} \item[a)] Verossimilhanças sem o efeito da covariável «»= sapply(Al.ml, logLik) @ \item[b)] Verossimilhanças com o efeito da covariável «»= sapply(Al.ml1, logLik) @ \end{itemize}

\noindent O maior valor obtido foi -45.16316 que corresponde ao modelo com a presença da covariável {\sf região} e para um valor de {\it kappa} de 3,5.

\subsubsection{Construção dos mapas temáticos}

Definindo o {\it grid} para a interpolação de valores por krigagem «results=hide»= gr ← pred_grid(Alum$borders, by=10) gr0← polygrid(gr, borders=Alum$borders, bound=T) ind.reg←numeric(nrow(gr0)) ind.reg[.geoR_inout(gr0,A1)]←1 ind.reg[.geoR_inout(gr0,A2)]←2 ind.reg[.geoR_inout(gr0,A3)]←3 ind.reg←as.factor(ind.reg) @

\noindent Efetuando a krigagem «results=hide»= set.seed(512) KC1←krige.control(trend.d=~regiao, trend.l=~ind.reg, obj.model=Al.ml1$Al.4) OC=output.control(n.pred=3) # define 3 simulacoes

pred←krige.conv(Alum, loc=gr, krige=KC1,output=OC) @ %$

\noindent A figura \ref{fig:mapas} é subdividida em quatro mapas. O primeiro trata-se de uma estimativa dos pontos com base na influência média de seus vizinhos, processo esse executado por uma krigagem convencional. Nos demais mapas os pontos são simulados com o modelo definido no ajuste de máxima-verossimilhança com seus respectivos parâmetros.

«fig=F, results=hide»= par(mfrow=c(2,2), mar=c(1,2,1,1)) image(pred, col=gray(seq(0.95,0.1,l=51))) image(pred, val=pred$simul[,1],col=gray(seq(0.95,0.1,l=51))) image(pred, val=pred$simul[,2],col=gray(seq(0.95,0.1,l=51))) image(pred, val=pred$simul[,3],col=gray(seq(0.95,0.1,l=51))) par(mfrow=c(1,1)) @ %$ \begin{figure}[!ht] \centering \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.99\textwidth} «fig=T, echo=F, results=hide»= par(mfrow=c(2,2), mar=c(1,2,1,1)) image(pred, col=gray(seq(0.95,0.1,l=51))) image(pred, val=pred$simul[,1],col=gray(seq(0.95,0.1,l=51))) image(pred, val=pred$simul[,2],col=gray(seq(0.95,0.1,l=51))) image(pred, val=pred$simul[,3],col=gray(seq(0.95,0.1,l=51))) par(mfrow=c(1,1)) @ %$ \vspace{-.5cm} \caption{Mapa de distribuição espacial da média (superior esquerdo) e de 3 simulações} \label{fig:mapas} \end{figure}