2

I have a problem with performing RF in R for a WorldView-2 Scene. I've never ever before used R for remote sensing classification, so I simply followed what is written here. I've loaded scene (8-band TIFF) with raster::brick and called it abcd, loaded simple 4-point, 2-class shapefile with OGR, run randomForest model and called it wv2rf and tried to use predict:

predict(abcd, wv2rf, "rftes.img", index=1, na.rm=TRUE)

R returned an error, that I cannot handle with:

Error in `[.data.frame`(blockvals, , f[j]) : undefined columns selected

I guessed that it has to do something with index value, but I tried to change it and nothing worked

whole code:

#loading data
abcd <- brick("13AUG20094646-M2AS-13EUSI-1283-01.tif")
sdata <- readOGR(dsn="mypath\\poligony testowe", layer="training")
v <- extract(x=abcd, y=sdata, df=TRUE)
sdata@data = data.frame(sdata@data, v[match(rownames(sdata@data), rownames(v)),]
sdata@data[3] <- NULL
#classification
wv2rf <- randomForest(x=sdata@data, y=(sdata@data[,"klasa"]), ntree=10, importance=TRUE)
predict(abcd, wv2rf, "rftes.img", index=1, na.rm=TRUE)

Here are files with my sdata, wv2rf and abcd variables, should work with load("filename") in R. Hope they help you help me.

http://speedy.sh/9f7cf/rf.zip

Community
  • 1
M'cin
  • 23
  • 5

1 Answers1

1

You are including your response variable in the training data. Besides not being a valid model this column is not available in your raster brick. Your "x" data needs to be indexed in the rf model so it excludes the response.

For example, if the first column is y ("klasa") and the rest are x: 

wv2rf <- randomForest(x=sdata@data[,2:ncol(sdata@data)], 
                      y=as.factor(sdata@data[,"klasa"]), 
                      importance=TRUE)

Based on your description, you only have 4 observations for discriminating 2 classes. Is this correct? If so, this is a quite invalid model! You need many more observations to apply a Bootstrap sampling model.

Jeffrey Evans
  • 31,750
  • 2
  • 47
  • 94
  • Yes, I just put 'anything' into the model to write script and check it, I will add valid training data later. Besides, it seems to work right now. Thanks a lot! – M'cin May 12 '14 at 22:15
  • One more question. Now I use polygons instead of points and I have problems with sdata@data = data.frame(sdata@data, match(rownames(sdata@data), rownames(v)),] sdata@data has only 26 rows, which represents number of polygons, while v has much more observations, that represents number of extracted pixels. Doing it this way leaves only 26 observations and names them incorrectly and I can't make it up on my own. – M'cin Aug 26 '14 at 19:10
  • This is really a different issue. Wen you extract points you get a value per observation. However, polygons result in a list object with all the pixel values intersecting each polygon. You can use a function like lapply() to return a statistic (e.g., mean, median, sd) that will correspond to the dimension of the sp polygon @data slot. – Jeffrey Evans Aug 26 '14 at 20:24