[GRASS-stats] Re: [GRASS-user] grass/stats mapping/prediction question

[Daniel McInerney]

Hi Andy,

 > I am unsure how to then move the 'outmap' back across to grass.
 > How can i convert the df to a spatial grid object.

AFAIK, 'predict' won't create a dataframe object.
R should return FALSE for is.data.frame(outmap) and
TRUE for is.numeric(outmap)

You can slot the model output to the AttributeList of
one of the SpatialGridDataFrames that you created when
you read in a GRASS raster and then use writeRAST6 to
write it back to GRASS.

e.g. using 'anmax' from your example

anmax$anmax <- outmap
#if that doesn't work, you might try
anmax$anmax <- as.numeric(outmap)
writeRAST6(anmax, "NameOfNewGRASSRaster", "anmax")

Regards,
Daniel.

cc: [hidden email]

andrew haywood wrote:

> Dear List,
>  
> I am having some problems analysing some ecoligical models in grass
> using the spgrass package through R.
>  
> I have 130 plot locations where i have observed presence/absence  of a
> species.
> I have followed a similar framework to the BUGSITE modelling
> example from Markus's 2003 grass gis handouts (Grass 5)
> I have no problems constructing the model based on the 130 plots and the
> environmental layers from grass.
>  
> However, I am having problems bringing all the maps through into R so I
> can make a prediction map.
> The region isnt too large 1600 by 800 cells at 10m resolution
> I can bring all the environmental layers through to R using readRAST6()
> which doesnt take too much time at all.
>  
> However i assume I must convert the spatial grid objects into dataframes
> to apply the predicted model function.
> So I then coerce them into dataframes using as.dataframe (this takes ages)
> I then merge all the dataframes into a single dataframe. (this takes ages)
>  
> I then apply the model predict to the new data frame.
>  
> I am unsure how to then move the 'outmap' back across to grass.
> How can i convert the df to a spatial grid object.
>  
> Im thinking i must be doing something wrong. As it quite quick to pull
> through the layers . But seems to take quite a lot of processing to get
> the layers into a datframe appropppriate for applying the predictions.
>  
> Any help would be greatly appreciated.
>  
> Andy
>  
>  
> #pull through environmental layers
> #FAST
> anmax <- readRAST6("anmax", ignore.stderr=TRUE)
> anmin <- readRAST6("anmin", ignore.stderr=TRUE)
> aspect <- readRAST6("aspect", ignore.stderr=TRUE)
> dem10_lidar <- readRAST6("dem10_lidar", ignore.stderr=TRUE)
>  
> #coerce to dataframe
> #SLOW
> mypred_anmaxDF<-as.data.frame(anmax)
> mypred_anminDF<-as.data.frame(anmin)
> mypred_aspectDF<-as.data.frame(aspect)
> mypred_dem10_lidarDF<-as.data.frame(dem10_lidar)
>  
> #merge into single dataframe
> #VERY SLOW
> merge_tmp<-merge(mypred_anmaxDF,mypred_anminDF)
> rm(mypred_anmaxDF,mypred_anminDF)
> merge_tmp1<-merge(merge_tmp,mypred_aspectDF)
> rm(merge_tmp,mypred_aspectDF)
> mypredDF<-merge(merge_tmp1,mypred_dem10_lidarDF)
>  
> #apply model
> outmap   <- predict(tree,newdata=mypredDF, type="class")
>  
>
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> grass-user mailing list
> [hidden email]
> http://lists.osgeo.org/mailman/listinfo/grass-user

_______________________________________________
grass-stats mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/grass-stats

[Roger Bivand]

On Wed, 19 Dec 2007, Daniel McInerney wrote:

The original questioner should have written to the grass-stats list in the
first place - thanks for CC-ing. See below for inline comments:

> Hi Andy,
>
>>  I am unsure how to then move the 'outmap' back across to grass.
>>  How can i convert the df to a spatial grid object.
>
> AFAIK, 'predict' won't create a dataframe object.
> R should return FALSE for is.data.frame(outmap) and
> TRUE for is.numeric(outmap)
>
> You can slot the model output to the AttributeList of

It hasn't been an AttributeList for a long time - the data slot *is* a
data frame. We *always* need the output of sessionInfo() to help -
update.packages() is most often very helpful too.

> one of the SpatialGridDataFrames that you created when
> you read in a GRASS raster and then use writeRAST6 to
> write it back to GRASS.
>
> e.g. using 'anmax' from your example
>
> anmax$anmax <- outmap
> #if that doesn't work, you might try
> anmax$anmax <- as.numeric(outmap)
> writeRAST6(anmax, "NameOfNewGRASSRaster", "anmax")
>
> Regards,
> Daniel.
>
> cc: [hidden email]
>
> andrew haywood wrote:
>>  Dear List,
>>
>>  I am having some problems analysing some ecoligical models in grass using
>>  the spgrass package through R.
>>
>>  I have 130 plot locations where i have observed presence/absence of a
>>  species. I have followed a similar framework to the BUGSITE modelling
>>  example from Markus's 2003 grass gis handouts (Grass 5) I have no
>>  problems constructing the model based on the 130 plots and the
>>  environmental layers from grass.

See the OSGeo tutorial September 2006:

http://www.foss4g2006.org/contributionDisplay.py?contribId=46&sessionId=59&confId=1

and the OSGeo Journal note:

http://www.osgeo.org/files/journal/final_pdfs/OSGeo_vol1_GRASS-R.pdf

for more up-to-date information.

>>
>>  However, I am having problems bringing all the maps through into R so I
>>  can make a prediction map.
>>  The region isnt too large 1600 by 800 cells at 10m resolution
>>  I can bring all the environmental layers through to R using readRAST6()
>>  which doesnt take too much time at all.
>>
>>  However i assume I must convert the spatial grid objects into dataframes
>>  to apply the predicted model function.

No, usually not at all, since the objects have a data.frame in the data
slot, and have the standard access methods.

>>  So I then coerce them into dataframes using as.dataframe (this takes ages)
>>  I then merge all the dataframes into a single dataframe. (this takes ages)
>>
>>  I then apply the model predict to the new data frame.
>>
>>  I am unsure how to then move the 'outmap' back across to grass.
>>  How can i convert the df to a spatial grid object.
>>
>>  Im thinking i must be doing something wrong. As it quite quick to pull
>>  through the layers . But seems to take quite a lot of processing to get
>>  the layers into a datframe appropppriate for applying the predictions.
>>
>>  Any help would be greatly appreciated.
>>
>>  Andy
>>
>>
>> # pull through environmental layers
>> # FAST
>>  anmax <- readRAST6("anmax", ignore.stderr=TRUE)
>>  anmin <- readRAST6("anmin", ignore.stderr=TRUE)
>>  aspect <- readRAST6("aspect", ignore.stderr=TRUE)
>>  dem10_lidar <- readRAST6("dem10_lidar", ignore.stderr=TRUE)

Wrong, do:

mydata <- readRAST6(c("anmax", "anmin", "aspect", "dem10_lidar"),
   ignore.stderr=TRUE)

then the data slot of the object is a data frame. Look at

summary(mydata)

for a sanity check.

>>
>> # coerce to dataframe
>> # SLOW
>>  mypred_anmaxDF<-as.data.frame(anmax)
>>  mypred_anminDF<-as.data.frame(anmin)
>>  mypred_aspectDF<-as.data.frame(aspect)
>>  mypred_dem10_lidarDF<-as.data.frame(dem10_lidar)
>>
>> # merge into single dataframe
>> # VERY SLOW
>>  merge_tmp<-merge(mypred_anmaxDF,mypred_anminDF)
>>  rm(mypred_anmaxDF,mypred_anminDF)
>>  merge_tmp1<-merge(merge_tmp,mypred_aspectDF)
>>  rm(merge_tmp,mypred_aspectDF)
>>  mypredDF<-merge(merge_tmp1,mypred_dem10_lidarDF)
>>
>>  #apply model

What is tree? You may need to do extra steps depending on what class(tree)
says - if you have used rpart() or some such, you may find that

outmap   <- predict(tree,newdata=mydata, type="class")

works,

or

outmap   <- predict(tree,newdata=as(mydata, "data.frame"), type="class")

Maybe just assign into mydata straight away:

mydata$outmap   <- predict(tree,newdata=as(mydata, "data.frame"),
   type="class")

given ?predict.rpart saying:

"If 'type="class"': (for a classification tree) a factor of
classifications based on the responses."

which looks like a vector for a vector response in the formula to rpart().
But do check what happens if there are NA in the newdata, because the
default predict() behaviour may be to drop those observations. Look at
summary(mydata).

Some formula-using model fitting functions just work, like lm() and
the predict() method for lm objects.

>From there, as Daniel wrote:

writeRAST6(mydata, "rpartpred", "outmap")

Hope this helps,

Roger

>>  outmap   <- predict(tree,newdata=mypredDF, type="class")
>>
>>
>>
>>  ------------------------------------------------------------------------
>>
>>  _______________________________________________
>>  grass-user mailing list
>>  [hidden email]
>>  http://lists.osgeo.org/mailman/listinfo/grass-user
> _______________________________________________
> grass-stats mailing list
> [hidden email]
> http://lists.osgeo.org/mailman/listinfo/grass-stats
>

--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: [phone]; fax [phone]
e-mail: [hidden email]
_______________________________________________
grass-stats mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/grass-stats

[Roger Bivand]

On Wed, 19 Dec 2007, Dylan Beaudette wrote:

> On Wednesday 19 December 2007, Roger Bivand wrote:
>> On Wed, 19 Dec 2007, Daniel McInerney wrote:
>>
>> The original questioner should have written to the grass-stats list in the
>>
>> first place - thanks for CC-ing. See below for inline comments:
>>> Hi Andy,
>>>
>>>>  I am unsure how to then move the 'outmap' back across to grass.
>>>>  How can i convert the df to a spatial grid object.
>>>
>>> AFAIK, 'predict' won't create a dataframe object.
>>> R should return FALSE for is.data.frame(outmap) and
>>> TRUE for is.numeric(outmap)
>>>
>>> You can slot the model output to the AttributeList of
>>
>> It hasn't been an AttributeList for a long time - the data slot *is* a
>> data frame. We *always* need the output of sessionInfo() to help -
>> update.packages() is most often very helpful too.
>>
>>> one of the SpatialGridDataFrames that you created when
>>> you read in a GRASS raster and then use writeRAST6 to
>>> write it back to GRASS.
>>>
>>> e.g. using 'anmax' from your example
>>>
>>> anmax$anmax <- outmap
>>> #if that doesn't work, you might try
>>> anmax$anmax <- as.numeric(outmap)
>>> writeRAST6(anmax, "NameOfNewGRASSRaster", "anmax")
>>>
>>> Regards,
>>> Daniel.
>>>
>>> cc: [hidden email]
>>>
>>> andrew haywood wrote:
>>>>  Dear List,
>>>>
>>>>  I am having some problems analysing some ecoligical models in grass
>>>> using the spgrass package through R.
>>>>
>>>>  I have 130 plot locations where i have observed presence/absence of a
>>>>  species. I have followed a similar framework to the BUGSITE modelling
>>>>  example from Markus's 2003 grass gis handouts (Grass 5) I have no
>>>>  problems constructing the model based on the 130 plots and the
>>>>  environmental layers from grass.
>>
>> See the OSGeo tutorial September 2006:
>>
>> http://www.foss4g2006.org/contributionDisplay.py?contribId=46&sessionId=59&
>> confId=1
>>
>> and the OSGeo Journal note:
>>
>> http://www.osgeo.org/files/journal/final_pdfs/OSGeo_vol1_GRASS-R.pdf
>>
>> for more up-to-date information.
>>
>>>>  However, I am having problems bringing all the maps through into R so I
>>>>  can make a prediction map.
>>>>  The region isnt too large 1600 by 800 cells at 10m resolution
>>>>  I can bring all the environmental layers through to R using readRAST6()
>>>>  which doesnt take too much time at all.
>>>>
>>>>  However i assume I must convert the spatial grid objects into
>>>> dataframes to apply the predicted model function.
>>
>> No, usually not at all, since the objects have a data.frame in the data
>> slot, and have the standard access methods.
>>
>>>>  So I then coerce them into dataframes using as.dataframe (this takes
>>>> ages) I then merge all the dataframes into a single dataframe. (this
>>>> takes ages)
>>>>
>>>>  I then apply the model predict to the new data frame.
>>>>
>>>>  I am unsure how to then move the 'outmap' back across to grass.
>>>>  How can i convert the df to a spatial grid object.
>>>>
>>>>  Im thinking i must be doing something wrong. As it quite quick to pull
>>>>  through the layers . But seems to take quite a lot of processing to get
>>>>  the layers into a datframe appropppriate for applying the predictions.
>>>>
>>>>  Any help would be greatly appreciated.
>>>>
>>>>  Andy
>>>>
>>>>
>>>> # pull through environmental layers
>>>> # FAST
>>>>  anmax <- readRAST6("anmax", ignore.stderr=TRUE)
>>>>  anmin <- readRAST6("anmin", ignore.stderr=TRUE)
>>>>  aspect <- readRAST6("aspect", ignore.stderr=TRUE)
>>>>  dem10_lidar <- readRAST6("dem10_lidar", ignore.stderr=TRUE)
>>
>> Wrong, do:
>>
>> mydata <- readRAST6(c("anmax", "anmin", "aspect", "dem10_lidar"),
>>    ignore.stderr=TRUE)
>>
>> then the data slot of the object is a data frame. Look at
>>
>> summary(mydata)
>>
>> for a sanity check.
>>
>>>> # coerce to dataframe
>>>> # SLOW
>>>>  mypred_anmaxDF<-as.data.frame(anmax)
>>>>  mypred_anminDF<-as.data.frame(anmin)
>>>>  mypred_aspectDF<-as.data.frame(aspect)
>>>>  mypred_dem10_lidarDF<-as.data.frame(dem10_lidar)
>>>>
>>>> # merge into single dataframe
>>>> # VERY SLOW
>>>>  merge_tmp<-merge(mypred_anmaxDF,mypred_anminDF)
>>>>  rm(mypred_anmaxDF,mypred_anminDF)
>>>>  merge_tmp1<-merge(merge_tmp,mypred_aspectDF)
>>>>  rm(merge_tmp,mypred_aspectDF)
>>>>  mypredDF<-merge(merge_tmp1,mypred_dem10_lidarDF)
>>>>
>>>>  #apply model
>>
>> What is tree? You may need to do extra steps depending on what class(tree)
>> says - if you have used rpart() or some such, you may find that
>>
>> outmap   <- predict(tree,newdata=mydata, type="class")
>>
>> works,
>>
>> or
>>
>> outmap   <- predict(tree,newdata=as(mydata, "data.frame"), type="class")
>>
>> Maybe just assign into mydata straight away:
>>
>> mydata$outmap   <- predict(tree,newdata=as(mydata, "data.frame"),
>>    type="class")
>>
>> given ?predict.rpart saying:
>>
>> "If 'type="class"': (for a classification tree) a factor of
>> classifications based on the responses."
>>
>> which looks like a vector for a vector response in the formula to rpart().
>> But do check what happens if there are NA in the newdata, because the
>> default predict() behaviour may be to drop those observations. Look at
>> summary(mydata).
>>
>> Some formula-using model fitting functions just work, like lm() and
>> the predict() method for lm objects.
>>
>>> From there, as Daniel wrote:
>>
>> writeRAST6(mydata, "rpartpred", "outmap")
>>
>> Hope this helps,
>>
>> Roger
>>
>>>>  outmap   <- predict(tree,newdata=mypredDF, type="class")
>>>>
>>>>
>
>
> An article on this in the OSGeo newsletter might be a nice way to document
> simple modeling examples with GRASS and R

Perhaps, and then there is the section in chapter 10 in "Open Source GIS
A GRASS GIS Approach", 3rd edition (my copy is still on its way, but from
the ToC, it looks as though pages 353-363 should be very helpful).

In fact, your site is a convenient collection of resources, I ought to
have mentioned it in my reply!

Roger

>
> D
>
>
>

--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: [phone]; fax [phone]
e-mail: [hidden email]
_______________________________________________
grass-stats mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/grass-stats