<html><head></head><body><div style="font-family: Verdana;font-size: 12.0px;"><div>
<div>
<div>Hello everyone,</div>
<div>For my Master thesis I use the GSIF Package for digital soil mapping by hands of a Random Forest model.<br/>
I have a <strong>spatial points data frame</strong> for the sampling locations and a <strong>spatial pixels data frame</strong> for the explanatory variables. I tried debugging and getting deeper into the error message but unfortunately I'm not able to solve the following problem.</div>
<div> </div>
<div><strong>Error in if (any(outside)) { : missing value where TRUE/FALSE needed</strong></div>
<div> </div>
<div>I’d be really thankful if somebody could have a glance at the following code. Perhaps the problem is really obvious….</div>
<div> </div>
<div>Thanks a lot</div>
<div> </div>
<div>Raphael</div>
</div>
<div> </div>
<div> </div>
<div>rm(list=ls())<br/>
library(spatial.tools)<br/>
library(nlme)<br/>
library(sp)<br/>
library(boot)<br/>
library(aqp)<br/>
library(plyr)<br/>
library(rpart)<br/>
library(splines)<br/>
library(gstat)<br/>
library(randomForest)<br/>
library(quantregForest)<br/>
library(plotKML)<br/>
library(GSIF)<br/>
library(rgdal)<br/>
library(sgeostat)<br/>
library(gstat)<br/>
library(rgeos)<br/>
library(maptools)<br/>
library(raster)<br/>
library(MASS)</div>
<div> </div>
<div> </div>
<div>#***********Get projection of the project********************</div>
<div><br/>
setwd("D:/STA/STA_R/pH_TA_20")<br/>
n <- "pH_Modelling"<br/>
suffix <- ".shp"<br/>
nn <- paste(n, suffix, sep = "")<br/>
projection <- readOGR(nn, n)<br/>
proj <- projection@proj4string</div>
<div>summary(proj)<br/>
#CRS arguments:<br/>
# +proj=somerc +lat_0=46.95240555555556 +lon_0=7.439583333333333 +k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel +units=m +no_defs</div>
<div><br/>
#******Import shapefile************</div>
<div><br/>
dataset <- readShapePoints(n,proj,verbose=TRUE,repair=FALSE)<br/>
summary(dataset)</div>
<div><br/>
# Object of class SpatialPointsDataFrame<br/>
# Coordinates:<br/>
# min max<br/>
# coords.x1 778617.174 778744.988<br/>
# coords.x2 205312.168 205428.537<br/>
# coords.x3 1527.428 1552.894<br/>
# Is projected: TRUE<br/>
# proj4string :<br/>
# [+proj=somerc +lat_0=46.95240555555556 +lon_0=7.439583333333333 +k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel +units=m +no_defs]<br/>
# Number of points: 60<br/>
# Data attributes:<br/>
# x y h Probe D0_5 D20_25 D50_55 D95_100 H_1 <br/>
# Min. :778617 Min. :205312 Min. :1527 g10S : 1 Min. :2.650 Min. :3.355 Min. :3.830 Min. :0.000 Min. :0.000 <br/>
# 1st Qu.:778646 1st Qu.:205343 1st Qu.:1536 g11S : 1 1st Qu.:3.169 1st Qu.:3.842 1st Qu.:4.242 1st Qu.:4.306 1st Qu.:0.000 <br/>
# Median :778674 Median :205367 Median :1541 g12S : 1 Median :3.375 Median :4.135 Median :4.513 Median :4.513 Median :2.720 <br/>
# Mean :778679 Mean :205368 Mean :1541 g13S : 1 Mean :3.945 Mean :4.301 Mean :4.602 Mean :4.569 Mean :1.586 <br/>
# 3rd Qu.:778715 3rd Qu.:205390 3rd Qu.:1548 g14S : 1 3rd Qu.:4.787 3rd Qu.:4.734 3rd Qu.:4.870 3rd Qu.:5.235 3rd Qu.:2.931 <br/>
# Max. :778745 Max. :205429 Max. :1553 g15S : 1 Max. :6.115 Max. :6.105 Max. :6.405 Max. :7.465 Max. :6.065 <br/>
# (Other):54 </div>
<div> </div>
<div>#**************Import Rasters of explaining variables****************</div>
<div><br/>
##Create raster Stack<br/>
AllRasters <- stack(list.files(path="D:/STA/STA_SAGA/DTMs/TAs_per_RasterSize/TA_20",pattern="*.tif$", full.names=TRUE))</div>
<div> </div>
<div>##Convert raster Stack to SPDF<br/>
Spatial_pixel_data_frame <- as(AllRasters, "SpatialPixelsDataFrame")</div>
<div> </div>
<div>##Set Projection of SPDF<br/>
projection(Spatial_pixel_data_frame) <- proj</div>
<div> </div>
<div>#*************Select some variables with no NAs********************<br/>
No_NAs <- Spatial_pixel_data_frame[c("dist2Flow_12c","dist2Flow_18c","dist2Flow_3c","dist2Flow_6c")]<br/>
summary(No_NAs)<br/>
# Object of class SpatialPixelsDataFrame<br/>
# Coordinates:<br/>
# min max<br/>
# x 778598.0 778752.2<br/>
# y 205295.4 205460.8<br/>
# Is projected: TRUE<br/>
# proj4string :<br/>
# [+proj=somerc +lat_0=46.95240555555556 +lon_0=7.439583333333333 +k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel +units=m +no_defs]<br/>
# Number of points: 637617<br/>
# Grid attributes:<br/>
# cellcentre.offset cellsize cells.dim<br/>
# s1 778598.1 0.2 771<br/>
# s2 205295.5 0.2 827<br/>
# Data attributes:<br/>
# dist2Flow_12c dist2Flow_18c dist2Flow_3c dist2Flow_6c <br/>
# Min. : 0.2843 Min. : 0.3726 Min. : 0.09251 Min. : 0.1514 <br/>
# 1st Qu.: 0.9141 1st Qu.: 0.9739 1st Qu.: 0.83619 1st Qu.: 0.8177 <br/>
# Median : 2.0820 Median : 2.0232 Median : 2.30100 Median : 2.2424 <br/>
# Mean :11.4551 Mean :11.4482 Mean :11.46038 Mean :11.4593 <br/>
# 3rd Qu.:16.2864 3rd Qu.:16.2419 3rd Qu.:16.30257 3rd Qu.:16.2985 <br/>
# Max. :83.1617 Max. :82.4677 Max. :84.16433 Max. :83.8265 </div>
<div> </div>
<div> </div>
<div>#************Run random forest Model in GSIF**************<br/>
omm <- fit.gstatModel(dataset, D20_25 ~ . , No_NAs , method="randomForest")</div>
<div> </div>
<div># Error in if (any(outside)) { : missing value where TRUE/FALSE needed<br/>
# In addition: Warning message:<br/>
# In .local(observations, formulaString, covariates, ...) :<br/>
# This method uses only 2D coordinates of the points. For 3D data consider using the 'geosamples-class'.</div>
<div> </div>
<div>options(error=recover)<br/>
</div>
<div> </div>
</div></div></body></html>