Não foi possível enviar o arquivo. Será algum problema com as permissões?
Diferenças

Diferenças

Aqui você vê as diferenças entre duas revisões dessa página.

Link para esta página de comparações

pessoais:exercicio_de_simulacao_geoestatistica_--_unioeste_13_11_07 [2008/07/15 17:34] (atual)
edson criada
Linha 1: Linha 1:
 +<​code>​ 
 +## ######################################​ 
 +## 
 +## SIMULACAO DE PROCESSO GEOESTATISTICO 
 +##      (Simulacao nao condicional) 
 +## 
 +##      Edson Antonio Alves da Silva 
 +##            (novembro, 2007) 
 +## ######################################​ 
 +##        Processo Gaussiano 
 +##    
 +##   ​Modelo:​ Y ~ N(mu ; Sigmasq R(phi)+Tausq I) 
 +## 
 +## ######################################​ 
 +## 
 +## Metodo manual usando os comando do R 
 +## 
 +## definindo os pontos amostrais 
 +grid <- expand.grid(X1=1:​3,​X2=1:​5) 
 +plot(grid) 
 +plot(grid, asp=1) 
 +class(grid) 
 +## 
 +## definindo a matriz de distâncias (euclidianas) 
 +mat.d <- round(as.matrix(dist(grid)),​2) 
 +mat.d 
 +## 
 +## definindo a matriz de covariancias 
 +## 
 +## ########################################​ 
 +## Fcao Correlacao: exponencial (-dist/​phi) 
 +phi <- 3 
 +sigma.sq <- 5 
 +mat.cov <- sigma.sq*(exp(-mat.d/​phi)) 
 +round(mat.cov,​2) 
 +## 
 +## Sorteando uma amostra normal multivariada com media mu 
 +## tomando k amostras univariadas e independentes N(0;1) 
 +## multiplicando pela matriz de covariancias e adicionando 
 +## a media (constante) 
 +mu <- 100 
 +set.seed(1956) 
 +Y <- rnorm(15)%*%mat.cov+mu 
 +range(Y) 
 +dados <- as.data.frame(cbind(grid,​Y=as.numeric(round(Y,​2)))) 
 +dados 
 +class(dados) 
 +##  
 +## convertendo os dados em um objeto geodata para uso da geoR 
 +## 
 +require(geoR) 
 +dados <- as.geodata(dados) 
 +class(dados) 
 +dados 
 +## ########################################​ 
 +## Fcao Correlacao: matern 
 +## 
 +kappa <- 0.5 
 +u.phi <- mat.d/phi 
 +mat.cov2 <- ifelse(mat.d > 0, 
 +                   ​(((2^(-(kappa-1)))/​gamma(kappa))*(u.phi^kappa)*besselK(x=u.phi,​nu=kappa)),​ 
 +                   1) 
 +set.seed(1956) 
 +Y2 <- rnorm(15)%*%mat.cov2+mu 
 +dados2 <- as.data.frame(cbind(grid,​Y=as.numeric(round(Y2,​2)))) 
 +dados2 
 +## 
 +## ###########################################################​ 
 +## 
 +## Usando a funcao GRF da GeoR 
 +## 
 +require(geoR) 
 +args(grf) 
 +## 
 +## Simulacao SEM tendencia (media constante) em grid REGULAR 
 +##  
 +## Definindo o grid (incrementado) 
 +## 
 +coord <- expand.grid(x1=(0:​9)*10+5,​x2=(0:​9)*10+5) 
 +plot(coord$x1,​coord$x2) 
 +A <- matrix(c(coord[,​1],​coord[,​2]),​ ncol=2, byrow=FALSE) 
 +pts1 <- expand.grid(seq(16,​24,​2),​seq(36,​44,​2)) 
 +pts2 <- expand.grid(seq(36,​44,​2),​seq(76,​84,​2)) 
 +pts3 <- expand.grid(seq(66,​74,​2),​seq(6,​14,​2)) 
 +pts4 <- expand.grid(seq(76,​84,​2),​seq(56,​64,​2)) 
 +B1 <- as.matrix(pts1) 
 +B2 <- as.matrix(pts2) 
 +B3 <- as.matrix(pts3) 
 +B4 <- as.matrix(pts4) 
 +area <- rbind(A,​B1,​B2,​B3,​B4) 
 +plot(area[,​1],​area[,​2],​pch=20,​cex=0.5) 
 +## 
 +## Parametros arbitrarios 
 +## 
 +## n=200 (dimensao do grid) 
 +## Fcao de correlacao esferica 
 +## Contribuicao = 15 
 +## Parametro de alcance phi = 60 
 +## Variacao de pequena escala tausq = 0 
 +## Parametro de trasformacao lambda = 1 (sem transformacao) 
 +## Distorcao de anisotropia:​ 
 +##       psi_A (angulo) = 0 
 +##       psi_B (razao entre eixos) = 1 
 +## 
 +dim(area) 
 +set.seed(1956) 
 +n=200 
 +sigma2 <- 15 
 +tausq <- 0 
 +phi <- 60 
 +Y3 <- grf(n,​grid=area,​cov.pars=c(sigma2,​phi),​ nugget=tausq,​cov.model='​sph'​) 
 +points(Y3,​col='​gray'​) 
 +## 
 +## 
 +## Simulacao SEM tendencia (media constante) em grid IRREGULAR 
 +## 
 +set.seed(1956) 
 +Y3.a <- grf(n,​grid='​irreg',​cov.pars=c(sigma2,​phi),​ nugget=tausq,​cov.model='​sph'​) 
 +plot(Y3.a$data,​Y3.a$coords[,​1]) 
 +points(Y3.a,​col='​gray'​) 
 +## 
 +## Simulacao COM tendencia em grid IRREGULAR 
 +## 
 +## Modelo: mu(x)= mu + beta1*coord(x1)+beta2*coord(x2) 
 +## 
 +## a) inicialmente geramos um modelo SEM tendencia 
 +set.seed(1956) 
 +Y3.b <- grf(n,​grid='​irreg',​cov.pars=c(sigma2,​phi),​ nugget=tausq,​cov.model='​sph'​) 
 +names(Y3.b) 
 +beta1 <- -5 
 +beta2 <- 5 
 +mu <- 10 
 +head(Y3.b$data) 
 +hist(Y3.b$data) 
 +set.seed(1956) 
 +Y3.b$data <- (Y3.b$data +mu)+beta1*Y3.b$coords[,​1]+beta2*Y3.b$coords[,​2] 
 +points(Y3.b) 
 +plot(Y3.b$data,​Y3.b$coords[,​2]) 
 +## 
 +## ##############################################################​ 
 +## 
 +## Processos nao gaussianos 
 +## 
 +## ##############################################################​ 
 +## 
 +## Distribuicao log-normal 
 +## 
 +Y4 <- grf(225,​grid='​reg',​ cov.pars=c(1,​0.24),​ lambda=0) 
 +hist(Y4$data) 
 +## 
 +## Distribuicao de Poisson 
 +## 
 +Y5 <- grf(225, grid="​reg",​ cov.pars=c(1,​ 0.25)) 
 +set.seed(1956) 
 +Y5$data <- rpois(225, exp(0.5 + Y5$data)) 
 +head(Y5$data) 
 +image(Y5, col=gray(seq(1,​0.2,​l=11))) 
 +text(Y5$coords[,​1],​ Y5$coords[,​2],​ Y5$data) 
 +barplot(table(Y5$data)) 
 +## 
 +## Distribuicao t-Student com 2 G.L. 
 +## 
 +## Simulacao SEM tendencia (media constante) em grid REGULAR 
 +## (nao foi utilizado a GRF da geoR) 
 +## Fcao Correlacao: exponencial (-dist/​phi) 
 +coord.t <- expand.grid(x1=(0:​9)*10+5,​x2=(0:​9)*10+5) 
 +dim(coord.t) 
 +mat.t <- round(as.matrix(dist(coord.t)),​2) 
 +phi <- 3 
 +sigma.sq <- 5 
 +mat.cov.t <- sigma.sq*(exp(-mat.t/​phi)) 
 +mu <- 100 
 +set.seed(1956) 
 +GL <- 1 
 +## O grau de liberdade (GL) define o decaimento da curva 
 +## Valores altos de GL equivalem a norma 
 +Y.t <- rt(100,​GL)%*%mat.cov.t+mu 
 +range(Y.t) 
 +hist(Y.t) 
 +dados.t <- as.data.frame(cbind(coord.t,​Y=as.numeric(Y.t,​2))) 
 +dados.t 
 +class(dados) 
 +</​code>​

QR Code
QR Code pessoais:exercicio_de_simulacao_geoestatistica_--_unioeste_13_11_07 (generated for current page)