[R-sig-Geo] how to export results of predict.randomForest (Edzer Pebesma)

Mon Dec 29 00:49:12 CET 2008

Hi Ann,
I've created prediction surfaces from random forest models and can try to help you. 

I don't know Rgdal well enough to guide you through your current approach. What you could do is create ascii grids for each of your raster/image bands, then call to those grids to create your prediction surface. 

Check out package = yaImpute and the AsciiGridPredict function within. I can send you some code to make this work if you'd like to try it. Just ask. Feel free to email off list, you can post the solution if it works.

Cheers,

Zack

> From: r-sig-geo-request at stat.math.ethz.ch
> Subject: R-sig-Geo Digest, Vol 64, Issue 25
> To: r-sig-geo at stat.math.ethz.ch
> Date: Sat, 27 Dec 2008 12:00:02 +0100
>
> Send R-sig-Geo mailing list submissions to
> r-sig-geo at stat.math.ethz.ch
>
> To subscribe or unsubscribe via the World Wide Web, visit
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> or, via email, send a message with subject or body 'help' to
> r-sig-geo-request at stat.math.ethz.ch
>
> You can reach the person managing the list at
> r-sig-geo-owner at stat.math.ethz.ch
>
> When replying, please edit your Subject line so it is more specific
> than "Re: Contents of R-sig-Geo digest..."
>
>
> Today's Topics:
>
> 1. Re: adding a background shape to ssplot (Edzer Pebesma)
> 2. Re: how to export results of predict.randomForest (Edzer Pebesma)
> 3. Re: adding a background shape to ssplot (Hengl, T.)
>
>
> ----------------------------------------------------------------------
>
> Message: 1
> Date: Fri, 26 Dec 2008 12:20:44 +0100
> From: Edzer Pebesma
> Subject: Re: [R-sig-Geo] adding a background shape to ssplot
> To: Jon Loehrke
> Cc: r-sig-geo at stat.math.ethz.ch
> Message-ID: 
> Content-Type: text/plain; charset=ISO-8859-1; format=flowed
>
> Jon Loehrke wrote:
>> Greetings and Happy Holidays!
>>
>>
>>
> ...
>> #However I can seem to pass the polygon colors argument through spplot
>> spplot(nc1, c("f","g"), col.regions=brewer.pal(10, "Set3"),
>> scales=list(draw = TRUE),
>> panel=function(...) {
>> panel.mapper()
>> sp.polygons(nc1)
>> }
>> )
>>
>>
> A better panel function would be this:
>
> spplot(nc1, c("f","g"), # col.regions=brewer.pal(10, "Set3"),
> scales=list(draw = TRUE),
> panel=function(...) {
> panel.mapper()
> panel.polygonsplot(...)
> }
> )
>
>
> see
>
> getMethod("spplot", "SpatialPolygonsDataFrame")
>
> to find out how the appropriate spplot method does it.
>
> I am worried, however, about the absence of the last 5 colors and factor
> levels f-j.
> --
> Edzer
>
>
>
> ------------------------------
>
> Message: 2
> Date: Fri, 26 Dec 2008 12:26:29 +0100
> From: Edzer Pebesma
> Subject: Re: [R-sig-Geo] how to export results of predict.randomForest
> To: Anne Axel
> Cc: r-sig-geo at stat.math.ethz.ch
> Message-ID: 
> Content-Type: text/plain; charset=ISO-8859-1; format=flowed
>
> I might have missed something obvious, but did you try something like
>
> bands$randomForest = output
> spplot(bands, "randomForest")
>
> ?
> --
> Edzer
>
>
> Anne Axel schrieb:
>> Hi Robert and list,
>>
>> Thanks so much for the Christmas message! I tried your great
>> suggestion and I got the following error:
>>
>>> spgrid <- SpatialGridDataFrame(bands, output)
>> Error in rbind(grid at cellcentre.offset, grid at cellcentre.offset +
>> (grid at cells.dim - :
>> no slot of name "cellcentre.offset" for this object of class
>> "SpatialGridDataFrame"
>>
>> So, I examined the structure of both files (see below). In addition,
>> there are no NAs in "imgtabla" or "bands", but there are NAs in the
>> rf.prediction ("output").
>>
>> I tried to writeGDAL(bands, "test2.tif") and this worked without any
>> errors. Could it have something to do with the fact that my output is
>> of class factor?
>>
>> Again, thanks Robert, for setting me on the right path. Any other
>> suggestions out there to clear this last hurdle?
>>
>> Best,
>> Anne Axel
>>
>>> str(bands)
>> Formal class 'SpatialGridDataFrame' [package "sp"] with 6 slots
>> ..@ data :'data.frame': 163982 obs. of 55 variables:
>> .. ..$ band1 : num [1:163982] 0 0 0 0 0 0 0 0 0 0 ...
>> .. ..
>> ..@ grid :Formal class 'GridTopology' [package "sp"] with 3 slots
>> .. .. ..@ cellcentre.offset: Named num [1:2] 0.5 0.5
>> .. .. .. ..- attr(*, "names")= chr [1:2] "x" "y"
>> .. .. ..@ cellsize : num [1:2] 1 1
>> .. .. ..@ cells.dim : int [1:2] 442 371
>> ..@ grid.index : int(0)
>> ..@ coords : num [1:2, 1:2] 0.5 441.5 0.5 370.5
>> .. ..- attr(*, "dimnames")=List of 2
>> .. .. ..$ : NULL
>> .. .. ..$ : chr [1:2] "x" "y"
>> ..@ bbox : num [1:2, 1:2] 0 0 442 371
>> .. ..- attr(*, "dimnames")=List of 2
>> .. .. ..$ : chr [1:2] "x" "y"
>> .. .. ..$ : chr [1:2] "min" "max"
>> ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots
>> .. .. ..@ projargs: chr " +proj=utm +zone=38 +south +a=6378137
>> +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
>>
>>> str(output)
>> Factor w/ 4 levels "dry","gallery",..: 3 3 3 3 3 3 3 3 3 3 ...
>>
>>
>>> class(bands)
>> [1] "SpatialGridDataFrame"
>> attr(,"package")
>> [1] "sp"
>>
>>
>>> class(output)
>> [1] "factor"
>>
>>> summary(output)
>> dry gallery notforest spiny NA's
>> 33009 5878 103772 21315 8
>>
>>
>>
>>
>>
>>
>>
>>
>>
>> At 12/24/2008 09:02 PM, you wrote:
>>> Hi Anne,
>>>
>>> The first argument of writeGDAL is "an object of class
>>> SpatialGridDataFrame-class or SpatialPixelsDataFrame-class"
>>> So you need to create a new SpatialGridDataFrame with the spatial
>>> properties of "multibands3.img" and with the RandomForest predictions
>>> ("output"). Note that there are some additional complications when
>>> there are NAs in your predictors (i.e. in "imgtabla"), but if there
>>> are none, I think that something like this might work:
>>>
>>> spgrid <- SpatialGridDataFrame(bands, output)
>>> writeGDAL(spgrid, "rf_results.tif")
>>>
>>> Robert
>>>
>>> On Thu, Dec 25, 2008 at 4:58 AM, Anne Axel
>>> wrote:
>>> Hi List,
>>>
>>> I am requesting help with exporting results from the
>>> predict.randomForest function. This is my first attempt at working
>>> with spatial data in R--I was able to import the data and run the
>>> randomForest model, but I cannot figure out how to export the data
>>> for use in ArcGIS or ENVI.
>>>
>>> I've pasted below the code I used for randomForest and exporting of
>>> results. If anyone has suggestions, I would be eternally grateful!
>>> Thanks, and Happy Holidays.
>>>
>>> Anne Axel
>>>
>>>
>>> Department of Fisheries and Wildlife
>>> Michigan State University
>>> 20 Natural Resources Building
>>> East Lansing, MI 48824 USA
>>> 517-930-2248
>>> axelanne at msu.edu
>>>
>>>
>>> ####Read in *img file with 55 bands
>>> bands=readGDAL("multibands3.img")
>>> imgtabla <- as(bands, "data.frame")
>>>
>>> ######Read in training data and run randomForest
>>> data=read.csv("bandnumsOutput_cond2.csv", header=TRUE)
>>> data2=data[,1:58]
>>> data.rf=randomForest(COVERTYPE~ ., data=data2, mtry=16, importance =
>>> TRUE, do.trace=100)
>>> print(data.rf)
>>> output=predict(data.rf, imgtabla)
>>>
>>> #### I've tried numerous ways of exporting the data, but I don't seem
>>> to have the data in a format that's acceptable to the functions. In
>>> addition, I've tried to coerce data into other formats, but nothing
>>> has worked. Any suggestions?
>>>
>>>
>>>> writeGDAL(bands, output, drivername = "GTiff", type = "Float32",
>>> mvFlag = NA, options=NULL)
>>> Error in .local(.Object, ...) :
>>> STRING_ELT() can only be applied to a 'character vector', not a
>>> 'integer'
>>> In addition: Warning messages:
>>> 1: In if (nchar(fname) == 0) stop("empty file name") :
>>> the condition has length> 1 and only the first element will be used
>>> 2: In if (nchar(filename) == 0) stop("empty file name") :
>>> the condition has length> 1 and only the first element will be used
>>>
>>>> writeGDAL(output, bands, drivername = "GTiff", type = "Float32",
>>> mvFlag = NA, options=NULL)
>>> Error in nchar(fname) : cannot coerce type S4 to character vector
>>>
>>>> write.ENVI(output)
>>> Error in dim(X) = c(nRow, nCol, nBand) : length-0 dimension vector is
>>> invalid
>>>
>>>> write.ascii.grid(output, "output", header = NULL, write.header = TRUE)
>>> Error: is.matrix(data) is not TRUE
>>>
>>>
>>>
>>>
>>>
>>> Anne C. Axel
>>> Department of Fisheries and Wildlife
>>> Michigan State University
>>> 20 Natural Resources Building
>>> East Lansing, MI 48824 USA
>>> 517-930-2248
>>> axelanne at msu.edu
>>>
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> R-sig-Geo at stat.math.ethz.ch
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>>
>> " Know your own bone. Gnaw at it, bury it, unearth it, and gnaw it
>> still."
>> -Henry David Thoreau
>>
>> "Outside of a dog, a book is a man's best friend. Inside of a dog,
>> it's too dark to read."
>> - Groucho Marx
>>
>> Anne C. Axel
>> Department of Fisheries and Wildlife
>> Michigan State University
>> 20 Natural Resources Building
>> East Lansing, MI 48824 USA
>> 517-930-2248
>> axelanne at msu.edu
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at stat.math.ethz.ch
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>
>
> ------------------------------
>
> Message: 3
> Date: Fri, 26 Dec 2008 15:01:27 +0100
> From: "Hengl, T."
> Subject: Re: [R-sig-Geo] adding a background shape to ssplot
> To: "Jon Loehrke" ,
> Message-ID: 
> Content-Type: text/plain
>
>
> take a look also at some examples from the book by Reimann et al. (2008). They actually load a jpg map (showing coast-line/borders) as a background plot (method "plotbg" in package "StatDA"):
>
> 6.4 Spatial Trends
> http://www.statistik.tuwien.ac.at/StatDA/R-scripts/page138.html
>
> 8.9 Data Subsets in Maps
> http://www.statistik.tuwien.ac.at/StatDA/R-scripts/page157.html
>
>
> HTH
>
> Tom Hengl
> http://spatial-analyst.net
>
> ________________________________
>
> From: r-sig-geo-bounces at stat.math.ethz.ch on behalf of Jon Loehrke
> Sent: sri 12/24/2008 7:01
> To: r-sig-geo at stat.math.ethz.ch
> Subject: [R-sig-Geo] adding a background shape to ssplot
>
>
>
> Greetings and Happy Holidays!
>
> I am an ~experienced R user but new to spatially oriented packages. I
> have been trying to add a 'background' shape to a ssplot. My example
> is to add a background state map to the sids data.
>
> My 2 questions, 1) What am I missing with the spplot - lattice
> interaction[this seems like it should be simple]?, and 2) is there a
> better way to add the 'background' map if the object is a
> SpatialPolygonsDataFrame?
>
> #Examples mostly from http://r-spatial.sourceforge.net/gallery/
> #Library's required for example
> library(sp)
> library(grid)
> library(lattice)
> library(maps)
> library(RColorBrewer)
> library(maptools)
>
> #First a Basic plot and an added map
>
> nc1 <- readShapePoly(system.file("shapes/sids.shp", package="maptools")
> [1], proj4string=CRS("+proj=longlat +datum=NAD27"))
>
> names(nc1)
> rrt <- nc1$SID74/nc1$BIR74
> brks <- quantile(rrt, seq(0,1,1/7))
> cols <- grey((length(brks):2)/length(brks))
> dens <- (2:length(brks))*3
> plot(nc1, col=cols[findInterval(rrt, brks, all.inside=TRUE)])
>
> # add US State map
> map('state', add=TRUE, xlim=par('usr')[1:2], ylim=par('usr')[3:4])
>
>
> #This works fine, but I have issues when using spplot.
>
> example
> nc1$f = factor(sample(1:5,100,replace=T),labels=letters[1:5])
> nc1$g = factor(sample(1:5,100,replace=T),labels=letters[6:10])
>
> spplot(nc1, c("f","g"), col.regions=brewer.pal(10, "Set3"),
> scales=list(draw = TRUE))
>
>
> # create a custom panel to add a map polygon to lattice plot
>
> panel.mapper <- function(...){
> library(maps)
> mp<-map('state', fill=TRUE, plot=FALSE)
> lpolygon(mp$x, mp$y, col='gray', ...)
> }
>
> #The panel mapper appears to work
> spplot(nc1, c("f","g"), col.regions=brewer.pal(10, "Set3"),
> scales=list(draw = TRUE), panel=function(...) {
> panel.mapper()
> })
>
> #However I can seem to pass the polygon colors argument through ssplot
> spplot(nc1, c("f","g"), col.regions=brewer.pal(10, "Set3"),
> scales=list(draw = TRUE),
> panel=function(...) {
> panel.mapper()
> sp.polygons(nc1)
> }
> )
>
>
> Thanks for the help and happiest holidays.
>
> Jon Loehrke
> Graduate Research Assistant
> Department of Fisheries Oceanography
> School for Marine Science and Technology
> University of Massachusetts
> 200 Mill Road, Suite 325
> Fairhaven, MA 02719
> jloehrke at umassd.edu
>> sessionInfo()
> R version 2.8.1 (2008-12-22)
> i386-apple-darwin8.11.1
>
> locale:
> en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8
>
> attached base packages:
> [1] grid stats graphics grDevices utils datasets
> methods base
>
> other attached packages:
> [1] lattice_0.17-17 RColorBrewer_1.0-2 maps_2.0-40
> maptools_0.7-15 foreign_0.8-29 sp_0.9-28
>
> loaded via a namespace (and not attached):
> [1] tools_2.8.1
>
>
> [[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
>
>
>
> [[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
>
>
> End of R-sig-Geo Digest, Vol 64, Issue 25
> *****************************************