<html><head><meta http-equiv=Content-Type content="text/html; charset=utf-8"><META name="Author" content="Novell GroupWise WebAccess"></head><body style='font-family: Helvetica, Arial, sans-serif; font-size: 13px; '>Jonathan,<br>Thanks again for your reply. <br>I'll reply inline ... please see below<br><br>>>>> Jonathan Greenberg <jgrn@illinois.edu> 03/26/14 10:55 AM >>><br>>Tim:<br>><br>>One thing to make sure is if the stack has the bands named properly<br>>>(they need to match the formula), e.g.:<br>>>names(mystack) <- c("b1","b2",...)<br><br>Yes, I think I'm good on layer names, as well as maintaining all the levels in categorical data (a classic gotcha for RF). <br><br>><br>>If you can also pass along the error you are getting, that will go a<br>>long way towards figuring out the issue! Cheers!<br><br>I was a little reluctant to pass on the error message until I knew a little more. The main reason being that the normal predict function with raster *does* work even though you might not expect it with the error message. Here is the error coming out of predict_rasterEngine:<br> <br>> sfQuickInit(cpus=6)<br>socket cluster with 6 nodes on host â€˜localhost’<br>> prediction_rf_prob <- predict_rasterEngine(object=rf.full, newdata=envStack, newtype='prob')<br>Error in predict.randomForest(object = object, newdata = newdata_df, mget(model_parameters)) : <br>  missing values in newdata<br>> sfQuickStop()<br> <br><br>Again, I can run these exact data using raster:::predict and I get good output.  To me that means I do not have 'missing values in newdata'. I welcome thoughts otherwise, however. The predict call is structured like this:<br><br>prediction_rf_prob <- predict(object=envStack, model=rf.full, type='prob', progress = 'text')<br><br>So, I've tried recreating the problem by tweaking the Tahoe data. Below, I've added a categorical layer to the brick. Unfortunately, it errors out with a new error. (Null external pointer).  (Am I right that raster:::predict and predict_rasterEngine use the term "object" differently?)<br><br>I'll paste a replacement set of code for the Tahoe data below our string.<br><br>I welcome comments on any of this!<br><br>Cheers, <br>Tim<br><br>>--j<br>><br>>On Wed, Mar 26, 2014 at 9:48 AM, Tim Howard <tghoward@gw.dec.state.ny.us> wrote:<br>>> Jonathan,<br>>> Thanks for your quick reply. I just tried your Tahoe example using x and y<br>>> in the call rather than a formula and it worked fine, so I was barking up<br>>> the wrong tree. Sorry. I'll try to delve deeper and put together a subset<br>>> example for you to check out if I can't get anywhere.<br>>><br>>> Thanks,<br>>> Tim<br>>><br>>>>>> Jonathan Greenberg <jgrn@illinois.edu> 3/26/2014 10:23 AM >>><br>>><br>>> Hi Tim:<br>>><br>>> Re: stack, yep, it should work. In general, you get a decreased<br>>> performance from stacks since you are having to read from multiple<br>>> files rather than a single one, but it will still benefit from<br>>> parallel processing.<br>>><br>>> Re: formula -- the best way for me to test this is to crop out a piece<br>>> of your image and send me the random forest model (use ?save), and the<br>>> call you were using. In theory you should be able to use anything you<br>>> would normally use on a data frame, but I'd have to play with it to<br>>> confirm! If you can set up those on e.g. google drive I can test it<br>>> out. Cheers!<br>>><br>>> --j<br>>><br>>> On Wed, Mar 26, 2014 at 6:49 AM, Tim Howard <tghoward@gw.dec.state.ny.us><br>>> wrote:<br>>>> Jonathan,<br>>>> Thank you for putting this together and for the example. I'm doing two<br>>>> things differently with randomForest ... I think perhaps one of them the<br>>>> function isn't handling.<br>>>><br>>>> First, based on recommendations from Andy Liaw (and ?randomForest), I<br>>>> don't<br>>>> use the formula interface but use x=<many columns>, y=<a column> in the<br>>>> call. Does predict_rasterEngine handle the absence of a formula in the<br>>>> object?<br>>>><br>>>> Second, I have many large rasters I want to run the predict on, so making<br>>>> a<br>>>> brick would be difficult. I use a rasterStack instead. Does your example<br>>>> work with a rasterStack?<br>>>><br>>>> I can dive deeper if any of this isn't clear or if these two tweaks work<br>>>> just fine for you. I was just trying to swap out this version of predict<br>>>> with another parallel version to evaluate speed and, while the alternate<br>>>> version works fine, predict_rasterEngine bailed on me.<br>>>><br>>>> Thanks in advance.<br>>>> Tim Howard<br>>>><br>>>><br>>>><br>>>>>>>>>><br>>>> Date: Tue, 18 Mar 2014 22:14:23 -0500<br>>>> From: Jonathan Greenberg <jgrn@illinois.edu><br>>>> To: "r-sig-geo@r-project.org" <R-sig-Geo@r-project.org><br>>>> Subject: [R-sig-Geo] Parallel predict now in spatial.tools<br>>>> Message-ID:<br>>>> <CABG0rfseg+p0h4HdYOK+_Za=OLMeTKHAT+TQn7g_FkEdYiunFQ@mail.gmail.com><br>>>> Content-Type: text/plain; charset=ISO-8859-1<br>>>><br>>>> R-sig-geo'ers:<br>>>><br>>>> I finally got around to building a parallel predict statement that<br>>>> I've included in version 1.3.7 (or later) of spatial.tools (check<br>>>> http://r-forge.r-project.org/R/?group_id=1492 for the status of the<br>>>> build), "predict_rasterEngine". It should, in theory, be a direct<br>>>> swap-in for the standard generic predict() statement. Currently, it<br>>>> will work on any predict.* statement that has the following features:<br>>>> 1) The data is passed to the predict as a data frame via a newdata<br>>>> parameter, and<br>>>> 2) The data is returned from the predict statement as a vector/matrix.<br>>>><br>>>> When using predict_rasterEngine, the object= parameter is your model,<br>>>> and the newdata= parameter is the raster/brick/stack to apply the<br>>>> model to on a pixel-by-pixel basis (note that the names of the layers<br>>>> must match the names of the predictor variables, in most cases).<br>>>><br>>>> I was hoping to get some stress-testing on this, since it is a fairly<br>>>> oft-requested function. If a predict.* function you'd like to use<br>>>> doesn't work, let me know which function it is (with some test data)<br>>>> and I'll see if I can tweak it to work.<br>>>><br>>>> Right now, I have confirmed this works with randomForest. Here's an<br>>>> example:<br>>>><br>>>> ######################<br>>>><br>>>> packages_required <- c("spatial.tools","doParallel","randomForest")<br>>>> lapply(packages_required, require, character.only=T)<br>>>><br>>>> # Load up a 3-band image:<br>>>> tahoe_highrez <- setMinMax(<br>>>> brick(system.file("external/tahoe_highrez.tif", package="spatial.tools")))<br>>>> tahoe_highrez<br>>>> plotRGB(tahoe_highrez)<br>>>><br>>>> # Load up some training points:<br>>>> tahoe_highrez_training_points <- readOGR(<br>>>> dsn=system.file("external", package="spatial.tools"),<br>>>> layer="tahoe_highrez_training_points")<br>>>><br>>>> # Extract data to train the randomForest model:<br>>>> tahoe_highrez_training_extract <- extract(<br>>>> tahoe_highrez,<br>>>> tahoe_highrez_training_points,<br>>>> df=TRUE)<br>>>><br>>>> # Fuse it back with the SPECIES info:<br>>>> tahoe_highrez_training_extract$SPECIES <-<br>>>> tahoe_highrez_training_points$SPECIES<br>>>><br>>>> # Note the names of the bands:<br>>>> names(tahoe_highrez_training_extract) # the extracted data<br>>>> names(tahoe_highrez) # the brick<br>>>><br>>>> # Generate a randomForest model:<br>>>> tahoe_rf <-<br>>>> randomForest(SPECIES~tahoe_highrez.1+tahoe_highrez.2+tahoe_highrez.3,<br>>>> data=tahoe_highrez_training_extract)<br>>>><br>>>> tahoe_rf<br>>>><br>>>> # This will run the predict in parallel:<br>>>> sfQuickInit()<br>>>> prediction_rf_class <-<br>>>><br>>>> predict_rasterEngine(object=tahoe_rf,newdata=tahoe_highrez,type="response")<br>>>> prediction_rf_prob <-<br>>>> predict_rasterEngine(object=tahoe_rf,newdata=tahoe_highrez,type="prob")<br>>>> sfQuickStop()<br>>>><br>>>> ###############<br>>>><br>>>> --j<br>>><br>>><br>>><br>>> --<br>>> Jonathan A. Greenberg, PhD<br>>> Assistant Professor<br>>> Global Environmental Analysis and Remote Sensing (GEARS) Laboratory<br>>> Department of Geography and Geographic Information Science<br>>> University of Illinois at Urbana-Champaign<br>>> 259 Computing Applications Building, MC-150<br>>> 605 East Springfield Avenue<br>>> Champaign, IL 61820-6371<br>>> Phone: 217-300-1924<br>>> http://www.geog.illinois.edu/~jgrn/<br>>> AIM: jgrn307, MSN: jgrn307@hotmail.com, Gchat: jgrn307, Skype: jgrn3007<br><br>#here is a script that adds a categorical layer to the tahoe brick, creates a RF model<br>#and then tries both raster:::predict and  predict_rasterEngine.<br><br>packages_required <- c("spatial.tools","doParallel","randomForest")<br>lapply(packages_required, require, character.only=T)<br><br># Load up a 3-band image:<br>tahoe_highrez <- setMinMax(<br>brick(system.file("external/tahoe_highrez.tif", package="spatial.tools")))<br>tahoe_highrez<br>plotRGB(tahoe_highrez)<br><br>#create a categorical layer from band 1<br>mat <- matrix(c(0,50,1,50,150,2,150,255,3),ncol=3,byrow=TRUE)<br>bnd1cat <- reclassify(tahoe_highrez[[1]], rcl=mat)<br>bnd1cat <- ratify(bnd1cat)<br>rat <- levels(bnd1cat)[[1]]<br>rat$types <- c('type1', 'type2', 'type3', 'type4')<br>rat$code <- c(1,2,3,4)<br>levels(bnd1cat) <- rat<br><br>#library(rasterVis) #if you want a categorical plot<br>#levelplot(bnd1_rc2)<br><br>tahoe_highrez <- addLayer(tahoe_highrez, bnd1cat)<br>names(tahoe_highrez) <- c("tahoeOne","tahoeTwo","tahoeThree","tahoeCatFour")<br><br># Load up some training points:<br>tahoe_highrez_training_points <- readOGR(<br>    dsn=system.file("external", package="spatial.tools"),<br>    layer="tahoe_highrez_training_points")<br><br># Extract data to train the randomForest model:<br>tahoe_highrez_training_extract <- extract(tahoe_highrez,tahoe_highrez_training_points,df=TRUE)<br><br># Fuse it back with the SPECIES info:<br>tahoe_highrez_training_extract$SPECIES <- tahoe_highrez_training_points$SPECIES<br><br># Note the names of the bands:<br>names(tahoe_highrez_training_extract) # the extracted data<br>names(tahoe_highrez) # the brick<br><br># convert to factor, ensure all the levels are there<br>tahoe_highrez_training_extract$tahoeCatFour <- factor(tahoe_highrez_training_extract$tahoeCatFour)<br>levels(tahoe_highrez_training_extract$tahoeCatFour) <- c(1,2,3,4)<br>str(tahoe_highrez_training_extract)<br><br># Generate a randomForest model:<br>tahoe_rf <- randomForest(y=tahoe_highrez_training_extract$SPECIES, <br>                            x=tahoe_highrez_training_extract[,2:5],<br>                            data=tahoe_highrez_training_extract)<br><br><br># try it with standard predict call -- this works                            <br>predict1_rf_prob <- predict(object=tahoe_highrez, model=tahoe_rf, type="prob")<br><br># try it with rasterEngine -- this fails<br>sfQuickInit()<br>predict2_rf_prob <- predict_rasterEngine(object=tahoe_rf,newdata=tahoe_highrez,type="prob")<br>sfQuickStop()<br><br>plot(predict1_rf_prob)<br><br><br><br></body></html>