# [R-sig-Geo] generate simulation data for a theoretical spatial model

Christoph Hofer christoph.hofer at env.ethz.ch
Thu Feb 4 08:55:33 CET 2010

rusers.sh,

You can simulate stationary and isotropic Gaussian random fields ( and transformation of it)  as well as max-stable random fields
by using the RandomFields package.
To get a conditional simulation of the correlated residuals of a linear model
you can first simulate an unconditional realization (\mu = 0 and a given covariance structure),
then choose the data location and make a simple kriging prediction of this realization.
Now, take the difference between the unconditional simulation and the simple kriging prediction and you
get one realization of a simulated conditional error term. Now, if you add this term to  a kriging estimation
you get one conditional simulation.

You find a detail descripitoin of this procedure called "Conditioning by kriging" in the book
Geostatistics: Modeling Spatial Uncertainty ( 1999) of Chilès, J.-P. and Delfiner, P

Sorry, for my bad english (I'm working on it).

Christoph

Am 03.02.2010 um 17:18 schrieb rusers.sh:

> For the method to generate the theoretical models, i agree with you
> completely. For binay (0/1) observations, it will be better to do as
> follows.
>  Y_i ~ B(p_i)
>  log(p_i/(1-p_i) = a1*x1+a2*x2+spatial effect
>  Take the above model as an example. Unfortunately, i still cannot figure
> out how to generate the conditional and unconditional spatial effects.
> #generate two covariates, x1 and x2
> library(mvtnorm)
> x1x2<-rmvnorm(100,mean=c(2,4),method="svd")
> colnames(x1x2)<-c("x1","x2")
> #Then, how to generate the conditional and unconditional spatial effects
> into the above model,respectively.
>  I think my mode of thinking on generating the spatial model in R must be
> at a dead end and needs guiding.
>  I appreciate your help very much.
>
>
> 2010/2/3 Paulo Justiniano Ribeiro Jr <paulojus at c3sl.ufpr.br>
>
>>
>> It is dificult if not irrealistic to set Y to be 0/1 (ou interger counts or
>> similar) in such model since this would impose severe contraints in the
>> a's and x's as well as in the model structure.
>>
>> This is why the hierarquical model structure is one possible
>> working around. The ideia is the same as in generalised linear models
>> relating the covariates (x's) and spatial effect to as function of the
>> expected value of Y instead of directly with Y.
>> In a "loose" notation:
>> Y_i ~ "some distribution" with E[Y_i] = \mu_i g(\mu_i) =
>> a1*x1+a2*x2+spatial effect
>>
>> where g() is a "convenient" function mapping (-Inf, +Inf)
>> to the parameter space of \mu_i
>>
>> Some examples:
>>
>> 1. For binay (0/1) observations a possible model would be
>> Y_i ~ B(p_i)
>> log(p_i/(1-p_i) = a1*x1+a2*x2+spatial effect
>>
>> 2. For count data:
>> Y_i ~ B(\lambda_i)
>> log(\lambda_i) = a1*x1+a2*x2+spatial effect
>>
>> 3. For Gaussian data
>> Y_i ~ N(\mu_i, \tau^2)
>> \mu_i = a1*x1+a2*x2+spatial effect
>> which in this particular case can be written as
>> Y_i = a1*x1+a2*x2+spatial effect
>>
>>
>>
>> Paulo Justiniano Ribeiro Jr
>> LEG (Laboratorio de Estatistica e Geoinformacao)
>> Caixa Postal 19.081
>> CEP 81.531-990
>> Curitiba, PR  -  Brasil
>> Tel: (+55) 41 3361 3573
>> Fax: (+55) 41 3361 3141
>> e-mail: paulojus AT  ufpr  br
>> http://www.leg.ufpr.br/~paulojus
>>
>>
>>
>> On Tue, 2 Feb 2010, rusers.sh wrote:
>>
>>  It works. The problem is that it only generates the simulated data based
>>> on our observed dataset,e.g. "meuse" here.
>>> I wonder if we can generate the simulated dataset from the user-specified
>>> model with covariates included, such as y~a1*x1+a2*x2+spatial effect. Y
>>> can
>>> be continuous or 0/1 variables. Something like this.
>>> The idea is we first specify a theoretical model, and then generate the
>>> simulated data based on this model. The coefficients and spatial effects
>>> are
>>> fixed by users, so we may study some new methods.
>>> Thanks.
>>>
>>> 2010/2/2 Edzer Pebesma <edzer.pebesma at uni-muenster.de>
>>>
>>>
>>>>
>>>> rusers.sh wrote:
>>>>
>>>> Hi Tomislav,
>>>>> Thanks for your info on unconditional simulation. For conditional
>>>>> simulations, i still cannot find any useful information.
>>>>> I searched the R site and didnot find the possible method to do
>>>>> conditional simulations.
>>>>> 1. CondSimu(RandomField): trend: Not programmed yet. (used by universal
>>>>> kriging)
>>>>> 2. grf(geoR): generates unconditional simulations of Gaussian random
>>>>> fields
>>>>> 3. sim.Krig(fields)  #Conditonal simulation of a spatial process
>>>>> It seems to be based on the actual dataset,not a theoretical model.
>>>>> 4. krige(gstat ):Simple, Ordinary or Universal, global or local, Point
>>>>> or
>>>>> Block Kriging,or simulation
>>>>> x <- krige(log(zinc)~x+y, meuse, meuse.grid, model = m, block =
>>>>> c(40,40),nsim=1)
>>>>>
>>>>>
>>>>
>>>> x <- krige(log(zinc)~x+y, meuse, meuse.grid, model = m, nmax=40, nsim=1)
>>>>
>>>> both adding the block=c(40,40) as well as omitting the nmax=40
>>>> tremendously
>>>> increased the computing time you needed, the second even more (in an
>>>> O(n^2)
>>>> manner) than the first.
>>>> --
>>>> Edzer
>>>>
>>>>
>>>>
>>>>
>>>> I used the above modified codes from krige(gstat ) example to see the
>>>>
>>>>> effect of "nsim", but unfortunately, it took a longer time and cannot
>>>>> get
>>>>> the results. I guess it used the simulation method to test the model,
>>>>> not
>>>>> what i want. (My system is XP, R2.10.0, gstat09.-64.)
>>>>> Anybody can give me further information on generating the conditional
>>>>> simulations from a theoretical model just like the unconditional
>>>>> examples
>>>>> that Tomislav provided?
>>>>> Thanks a lot.
>>>>>
>>>>>
>>>>> 2010/1/31 Tomislav Hengl <hengl at spatial-analyst.net>
>>>>>
>>>>>
>>>>>
>>>>> Dear rusers.sh,
>>>>>>
>>>>>> Here are few simple examples of how to simulate (not-normal)
>>>>>> distributions and point processes using geoR and spatstat:
>>>>>>
>>>>>> http://spatial-analyst.net/book/node/388
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> http://leg.ufpr.br/geoR/geoRdoc/vignette/geoRintro/geoRintrose8.html#x9-120008
>>>>>>
>>>>>> I guess that covariates can be also included (I guess that you then
>>>>>> need
>>>>>> to switch to conditional simulations - not sure).
>>>>>>
>>>>>> This should also work for lattice (polygon) data so that you will have
>>>>>> jumps in values (but I guess you would still work in gridded systems?).
>>>>>>
>>>>>> T. Hengl
>>>>>> http://home.medewerker.uva.nl/t.hengl/
>>>>>>
>>>>>>
>>>>>> rusers.sh wrote:
>>>>>>
>>>>>>
>>>>>>
>>>>>> Hi all,
>>>>>>> In classical statistics, we always need to generate a theoretical
>>>>>>> model
>>>>>>> such as y=a+b1*x1+b2*x2+e to study some new estimation content. I am
>>>>>>> wondering how to generate the similar spatial dataset for a
>>>>>>> theoretical
>>>>>>> model.
>>>>>>> Say y is response variable, x1 and x2 are explanatory variables.
>>>>>>> 1. If y is a continous variable, how should we generate the dataset
>>>>>>> for
>>>>>>> a
>>>>>>> theoretical spatial point process model in R?
>>>>>>> 2. If y is a continous variable, how should we generate the dataset
>>>>>>> for
>>>>>>> a
>>>>>>> theoretical spatial lattice data model in R?
>>>>>>> 3. If y is 0/1 binary variable, how should we generate the dataset for
>>>>>>> a
>>>>>>> theoretical spatial point process model in R?
>>>>>>> 4. If y is 0/1 binary variable, how should we generate the dataset for
>>>>>>> a
>>>>>>> ttheoretical spatial lattice data in R?
>>>>>>> spatstat and other packages allow us to generate a dataset of a
>>>>>>> specified
>>>>>>> point process and other models, but it seems that they donot allow us
>>>>>>> to
>>>>>>> include possible explanatory variables into a theoretical model. Maybe
>>>>>>> i
>>>>>>> missed some ideas in them.
>>>>>>> Anybody can express some ideas or point out some useful resources on
>>>>>>> the
>>>>>>> above four different situations? Small examples in R are preferred.
>>>>>>> Thanks a lot.
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> _______________________________________________
>>>>>> R-sig-Geo mailing list
>>>>>> R-sig-Geo at stat.math.ethz.ch
>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>> --
>>>> Edzer Pebesma
>>>> Institute for Geoinformatics (ifgi), University of M?nster Weseler Stra?e
>>>> 253, 48151 M?nster, Germany. Phone: +49 251 8333081, Fax: +49 251 8339763
>>>>
>>>> http://ifgi.uni-muenster.de http://www.52north.org/geostatistics
>>>> e.pebesma at wwu.de
>>>>
>>>>
>>>>
>>>
>>> --
>>> -----------------
>>> Jane Chang
>>> Queen's
>>>
>>>       [[alternative HTML version deleted]]
>>>
>>>
>>>
>
>
> --
> -----------------
> Jane Chang
> Queen's
>
> 	[[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo