[R-sig-Geo] imaging geoRglm binomial krig

Ken Nussear knussear at mac.com
Thu Mar 5 13:49:35 CET 2009


Hi Noicola,
Thanks for the reply

Mine were originally made in ArcMap by a colleague of mine, then  
imported into Grass, which I use. The point locations (x,y) and the  
covariates are all in a vector attribute table for the points and all  
values are pulled direclty from a sqlite backend via RSQLite (Code  
Below). I haven't tried pred_grid. Ill give it a shot.

drv <- dbDriver("SQLite")
con <- dbConnect(drv, dbname="/Users/knussear/Documents/Research/GIS/ 
grassdata/DTSpatial/Working/sqlite.db")
dtdata2 <- dbSendQuery(con, "SELECT * from LDS1kGridE WHERE XMIN >=  
415311.9 and XMIN <= 555311.9 and YMIN >= 3823651 and YMIN <= 3918651")
kriglocs <- fetch(dtdata2, n=-1)

#AnUnits='WM' and
kriglocs $TRAN_LNTH <- kriglocs $TRAN_LNTH/1000
kriglocs $MDEP <- kriglocs $MDEP/1000
kriglocs $HWYS_Dist2 <- kriglocs $HWYS_Dist2/1000
kriglocs $pop <- kriglocs $Pop_15km/1000
kriglocs $ER <- kriglocs $LIVE_COUNT/kriglocs $TRAN_LNTH

WMgeo.kriglocs <- as.geodata(kriglocs, coords.col=c(5,7), data.col=44,  
covar.col=c(34,37,42,43))


Ken


On Mar 5, 2009, at 3:18 AM, Nicola Batchelor wrote:

> Hi Ken,
>
> I've also been using geoRglm for similar analyses, and found the  
> prediction
> to be a bit fickle when it comes to the prediction grids!
>
> I was having similar problems to yours and the only thing that I  
> could do to
> make it work was to alter the way I made my prediction grids.  How  
> are you
> creating yours?
>
> What I was doing initially was creating the prediction points in  
> ArcMap,
> adding the covariate data to it and then exporting and importing to  
> R which
> didn't work.
>
> I ended up changing my method, so I'd create the prediction points  
> in R
> using the pred_grid function, then export that and import the points  
> to
> ArcMap.  Then I'd attach all the covariate data to the prediction  
> points and
> then export from ArcMap and import to R as a covariate object (and  
> also use
> the prediction points created using pred_grid as my prediction  
> locations).
>
> I hope this helps you a bit, but perhaps your problem is not the  
> same as
> mine was!
>
> Nicola
>
>
> Date: Wed, 04 Mar 2009 11:07:16 -0800
> From: Ken Nussear <knussear at mac.com>
> Subject: [R-sig-Geo] imaging geoRglm binomial krig
> To: r-sig-geo at stat.math.ethz.ch
> Message-ID: <1898D59A-04AF-4625-B19C-E8DACA88A133 at mac.com>
> Content-Type: text/plain; charset=US-ASCII; format=flowed; delsp=yes
>
> Hi
>
> Wondering if anyone has had success imaging a krige produced  using
> geoRglm?
>
> The dataset has 3467 locations....
>
> Here is the code used to build the krige and attempted image.
>
> itdsglm1 <- trend.spatial(~TRAN_LNTH + Maxent + HWYS_Dist2 + MDEP +
> pop, WMgeo.Sign)
>
> itlsglm1 <- trend.spatial(~TRAN_LNTH + Maxent + HWYS_Dist2 + MDEP +
> pop, WMgeo.kriglocs)
>
> kmod <- model.glm.control(trend.d = itdsglm1, trend.l = itlsglm1,
> cov.model = "exponential")
>
> kprior <- prior.glm.control(beta.prior = "flat", sigmasq.prior =
> "fixed", sigmasq= 0.0307, phi.prior = "fixed", phi = 7748,
> phi.discrete = NULL, tausq.rel = 5.276)
>
>
> bkcontrl <- mcmc.control(S.scale = .07, thin=40, phi.start= 7748,
> phi.scale = .2, burn.in=10)
> kgout <- output.glm.control(sim.posterior=T, sim.predict=T,
> keep.mcmc.sim=T, inference=T, messages=T, quantile=c(0.25,0.5,0.975))
>
>
> bkb <- binom.krige.bayes(WMgeo.Sign, locations = WMgeo.kriglocs
> $coords,  model= kmod, mcmc.input= bkcontrl, prior=kprior)
>
> image(bkb, locations=WMgeo.kriglocs$coords,
> values.to.plot='simulation', number.col=1, messates=T)
>
>
>
> image.glm.krige.bayes(bkb, locations=WMgeo.kriglocs$coords,
> values.to.plot=c("median"), number.col=1, coords.data=WMgeo.kriglocs
> $coords, x.leg=c(410312,561312), y.leg=c(3822651,3920651), messages=T)
>
>
> I get the following errors and I can't figure out what it wants....
>
>
> mapping the medians of the predictive distribution
> Error in image.default(x = c(410311.919949, 411311.919949,
> 412311.919949,  :
>  dimensions of z are not length(x)(-1) times length(y)(-1)
> In addition: Warning messages:
> 1: In if (coords.data) coords.data <- eval(attr(x,  
> "data.locations")) :
>  the condition has length > 1 and only the first element will be used
> 2: In matrix(values.loc, ncol = ny) :
>  data length [3467] is not a sub-multiple or multiple of the number
> of rows [36]
>
>
>
> Ken
>
>
>
> -- 
> The University of Edinburgh is a charitable body, registered in
> Scotland, with registration number SC005336.
>



More information about the R-sig-Geo mailing list