[R-sig-Geo] Using GStat for kriging with Anisotropy

[Dan Turenne]

Hello,


I had considered this possibility but I was concerned about the distortion of distances.  If the coordinates are being changed from a projection on a sphere to a 2 dimensional Cartesian plane won't the distances between points be affected?  The area that I am interpolating across is fairly large and has a diagonal of roughly 850 km.  I apologize if these are basic questions, I am writing my thesis for my master's degree in actuarial science  and this project is my first exposure to the field of geostatistics.


Many thanks,

Daniel Turenne



________________________________
From: R-sig-Geo <r-sig-geo-bounces at r-project.org> on behalf of Edzer Pebesma <edzer.pebesma at uni-muenster.de>
Sent: June 15, 2016 3:49 PM
To: r-sig-geo at r-project.org
Subject: Re: [R-sig-Geo] Using GStat for kriging with Anisotropy

Dan, have you considered projecting your data to Carthesian coordinates
(using spTransform), doing the modelling and interpolation there, and
then back-transform results to long/lat?

On 15/06/16 20:00, Dan Turenne wrote:
> Hello,
>
>
> I am trying to use kriging to make temperature predictions but I have run into a slight issue.  I would like to include anisotropy in my variogram model but GStat does not support anisotropy when coordinates are given in longitude and latitude.  Here is a example of what my code looks like, where dailyData is a data frame containing long/lat coordinates and the mean temperature data:
>
>
>
>     coordinates(dailyData) = ~long+lat
>
>     proj4string(sp) = CRS("+init=epsg:4326")
>
>
>     anisParams=estimateAnisotropy(dailyData,"meanTemp")
>
>     sample.variogram=variogram(meanTemp~1,dailyData)
>
>
>     sillEst=var(dailyData$meanTemp)
>
>     nuggetEst=0
>
>     rangeEst=200
>
>
>     anis.variogram=vgm(nugget=nuggetEst, model="Exp", range=rangeEst, psill=sillEst,
>
>                                           anis=c(anisParams$direction,1/anisParams$ratio))
>
>     anis.fitVariogram=fit.variogram(sample.variogram,anis.variogram)
>
>
>     surface=krige(meanTemp~1,locations=dailyData,newdata=fData,model=anis.fitVariogram)
>
>
>
> However running this code produces the following error:
>
>
>     Error in predict.gstat(g, newdata = newdata, block = block, nsim = nsim,  :
>     value not allowed for: nnot be defined for long/lat data, anisotropy cannot be defined
>
>
> I know that the problem is that GStat can't model anisotropy while using great circle distances, but I am wondering if there is a workaround that exists which allows anisotropy to be modelled while using long/lat coordinates?  All of my long/lat data is in decimal format (i.e lat=45.62, long=-75.48).  I have searched through the R-Sig-Geo mailing list and have seen this issue mentioned but I have not been able to find a solution that works for my code.  Any help on this subject would be greatly appreciated.
>
>
> Thank you in advance,
>
>
> Daniel Turenne
>
> University of Manitoba
>
>        [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
R-sig-Geo - R Special Interest Group on using Geographical ...<https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
stat.ethz.ch
R-sig-Geo -- R Special Interest Group on using Geographical data and Mapping About R-sig-Geo



>

--
Edzer Pebesma
Institute for Geoinformatics  (ifgi),  University of Münster
Heisenbergstraße 2, 48149 Münster, Germany; +49 251 83 33081
Journal of Statistical Software:   http://www.jstatsoft.org/
Computers & Geosciences:   http://elsevier.com/locate/cageo/
Spatial Statistics Society http://www.spatialstatistics.info


	[[alternative HTML version deleted]]