Diferenças

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

--- artigos:ernesto3:simpj [2008/10/09 13:26] – ernesto
+++ artigos:ernesto3:simpj [2008/10/16 00:50] (atual) – paulojus
@@ Linha 1: / Linha 1: @@
 ====== A joint model proposal???? ======
+  - {{artigos:ernesto3:joint.rnw|Arquivo Sweave com o modelo proposto}}
 ===== Ideias for simulation and model =====
@@ Linha 76: / Linha 77: @@
 ===== Alternativa usando CDA =====
-</code>
+<code R>
-gs <- expand.grid((0:10)/10, (0:10)/10)
+require(lattice)
+require(MASS)
+require(geoR)
+gs <- expand.grid((0:10)/10, (0:10)/10)
 # basic gaussians to be used posteriorly
 s2 <- 0.5
@@ Linha 85: / Linha 90: @@
 Sig[4,1] <- 0.9
 Sig <- Sig*s2
+Sig
 m0 <- mvrnorm(nrow(gs), c(2,0,0,0,0), Sig)
+plot(as.data.frame(m0))
+cor(m0)
-# Generate a spatial gaussian process Y
+# Generate a spatial Gaussian process Y
 phi <- 0.25
 sigmasq <- 1
-Y <- grf(100, grid=gs, cov.pars=c(sigmasq, phi))
+Y <- grf(grid=gs, cov.pars=c(sigmasq, phi))
 Y$data <- exp(m0[,1]+Y$data)
+var(log(Y$data))
 ## option 1: age compositions independent from Y
 # Generate age compositions independent from Y
 comp.1 <- t(apply(cbind(exp(m0[,c(2,3)]),1),1,function(x) x/sum(x)))
+dim(comp.1)
+apply(comp.1, 2, range)
+cor(cbind(comp.1, log(Y$data)))
+plot(as.data.frame(cbind(comp.1, log(Y$data))))
+lr.1 <- data.frame(lr13 = log(comp.1[,1]/comp.1[,3]), lr23 = log(comp.1[,2]/comp.1[,3]), Y=log(Y$data))
+cor(lr.1)
+cor(lr.1, met="spea")
+plot(lr.1)
 # Build abundance at age
 Yi.1 <- comp.1*Y$data
+dim(Yi.1)
+cor(cbind(Yi.1, log(Y$data)))
+plot(as.data.frame(cbind(Yi.1, log(Y$data))))
 ## option 2: age compositions dependent from Y
 # Generate age compositions dependent from Y
 comp.2 <- t(apply(cbind(exp(m0[,c(3,4)]),1),1,function(x) x/sum(x)))
+apply(comp.1, 2, range)
+cor(cbind(comp.2, log(Y$data)))
+cor(cbind(comp.2, log(Y$data)), met="spea")
+plot(as.data.frame(cbind(comp.2, log(Y$data))))
-# Build abundance at age
+lr.2 <- data.frame(lr13 = log(comp.2[,1]/comp.2[,3]), lr23 = log(comp.2[,2]/comp.2[,3]), Y=log(Y$data))
-Yi.2 <- comp.2*Y$data
+cor(lr.2)
+cor(lr.2, met="spea")
+plot(lr.2)
-df0 <- data.frame(abund=c(Yi.1,Yi.2), comp=c(comp.1,comp.2),Y=Y$data, opt=rep(c("indep","dep"),rep(length(comp.1),2)), age=rep(rep(c(1,2,3),2),rep(length(comp.1)/3,6)))
+lr.2a <- data.frame(lr12 = log(comp.2[,1]/comp.2[,2]), lr32 = log(comp.2[,3]/comp.2[,2]), Y=log(Y$data))
+cor(lr.2a)
+cor(lr.2a, met="spea")
+plot(lr.2a)
+# Build abundance at age
+Yi.2 <- comp.2*Y$data
+dim(Yi.2)
+cor(cbind(Yi.2, log(Y$data)))
+plot(as.data.frame(cbind(Yi.2, log(Y$data))))
+df0 <- data.frame(abund=c(Yi.1,Yi.2), comp=c(comp.1,comp.2),Y=Y$data, opt=rep(c("indep","dep"),
+                  rep(length(comp.1),2)), age=rep(rep(c(1,2,3),2),rep(length(comp.1)/3,6)))
 # plots
 xyplot(comp~log(Y)|age*opt,data=df0)
@@ Linha 116: / Linha 156: @@
 O problema que temos nestes casos em que simulamos Y e depois obtemos Yi por Y*pi é que a estrutura espacial dos Yi vai ser exactamente a mesma. Se acrescentarmos variabilidade com estrutura espacial temos que recalcular os pi ...
+PJ: De fato mas para fazer diferente teria que simular de outra forma pois aqui só tem um Y geoestatistico gerado. O modelo multivariado  conjunto para os LR seria a alternativa. Por outro lado eu nao vejo muito problema nesta fato pois as correlacoes entra as composicoes de certa forma tratam isto.

Diferenças

Navegação

Busca

Ferramentas

QR Code