You have to dive a little bit into what factor() does. One or more of your classes never has the highest probability, it seems:
> factor(c(1,2,3,5)) [1] 1 2 3 5 Levels: 1 2 3 5 > factor(c(1,2,3,5), labels = c("1", "2", "3", "4", "5")) Error in factor(c(1, 2, 3, 5), labels = c("1", "2", "3", "4", "5")) : invalid labels; length 5 should be 1 or 4 > levels(factor(c(1,2,3,5))) [1] "1" "2" "3" "5" > factor(c(1,2,3,5), levels=1:5, labels = c("1", "2", "3", "4", "5")) [1] 1 2 3 5 Levels: 1 2 3 4 5 > factor(c(1,2,3,5), levels=1:5, labels = c("a", "b", "c", "d", "e")) [1] a b c e Levels: a b c d e On 07/24/2011 08:52 PM, Matevž Pavlič wrote: > Hi Edzer, > > thanks for the help. It works ok on the sampleI provided, even if I enlarge > the number of classes to 10. But when i try to use it on my dataset which has > identical names of predictions and number of predisctions and covarinaces I > get an error : > > Error in factor(apply(cbind(z$nic.pred, z$ena.pred, z$dva.pred, z$tri.pred, > : > invalid labels; length 10 should be 1 or 8 > > This is the oputput of the prediction whit > >> str(z): > Formal class 'SpatialPixelsDataFrame' [package "sp"] with 7 slots > ..@ data :'data.frame': 8736 obs. of 65 variables: > .. ..$ devet.pred : num [1:8736] 0.001153 0.000876 0.100325 0.101093 > 0.102047 ... > .. ..$ devet.var : num [1:8736] 0.00548 0.00548 0.00549 0.0055 0.00551 > ... > .. ..$ dva.pred : num [1:8736] 0.182 0.185 0.177 0.194 0.191 ... > .. ..$ dva.var : num [1:8736] 0.0744 0.0744 0.0745 0.0752 0.0753 ... > .. ..$ ena.pred : num [1:8736] 0.2009 0.1895 0.0476 0.0522 0.0482 ... > .. ..$ ena.var : num [1:8736] 0.149 0.149 0.15 0.155 0.155 ... > .. ..$ nic.pred : num [1:8736] 0.139 0.142 0.206 0.0252 0.0234 ... > .. ..$ nic.var : num [1:8736] 0.258 0.259 0.26 0.267 0.268 ... > .. ..$ osem.pred : num [1:8736] 9.15e-05 0.00 5.56e-05 6.82e-04 > 1.13e-03 ... > .. ..$ osem.var : num [1:8736] 0.00381 0.00382 0.00382 0.00386 > 0.00386 ... > .. ..$ pet.pred : num [1:8736] 0.000115 0.000303 0.000875 0.000423 > 0.000345 ... > .. ..$ pet.var : num [1:8736] 0.0126 0.0126 0.0126 0.0133 0.0133 ... > .. ..$ sedem.pred : num [1:8736] 0.0934 0.0924 0.0756 0.1084 0.1157 ... > .. ..$ sedem.var : num [1:8736] 0.0848 0.085 0.0852 0.0874 0.0877 ... > .. ..$ set.pred : num [1:8736] 0.271 0.278 0.27 0.395 0.389 ... > .. ..$ set.var : num [1:8736] 0.0881 0.0882 0.0885 0.0902 0.0904 ... > .. ..$ stiri.pred : num [1:8736] 0.001931 0.000777 0.0046 0 0 ... > .. ..$ stiri.var : num [1:8736] 0.101 0.101 0.101 0.106 0.106 ... > .. ..$ tri.pred : num [1:8736] 0.11 0.111 0.118 0.123 0.129 ... > .. ..$ tri.var : num [1:8736] 0.134 0.134 0.134 0.136 0.136 ... > > I would like to create a map from predictions of course. > > What i found out is that your latest code works until labels line gets > called.But funny thing is , when I use this line : > z$class = factor(apply(cbind(z$nic.pred,z$ena.pred, z$dva.pred, z$tri.pred, > z$stiri.pred, z$pet.pred, z$set.pred, z$sedem.pred, z$osem.pred, > z$devet.pred), 1, which.max)) > > i get : > >> levels(z$class) > [1] "1" "2" "3" "4" "5" "6" "7" "8" > > Any ideas? > > tnx, m > -----Original Message----- > From: r-sig-geo-boun...@r-project.org > [mailto:r-sig-geo-boun...@r-project.org] On Behalf Of Edzer Pebesma > Sent: Sunday, July 24, 2011 8:04 PM > To: r-sig-geo@r-project.org > Subject: Re: [R-sig-Geo] Indicator kriging map creation > > Matevz, try > > z$class = factor(apply(cbind(z$dva.pred, z$ena.pred, z$nic.pred, > z$stiri.pred, z$tri.pred), 1, which.max), > labels=c("dva", "ena", "nic", "stiri", "tri")) > spplot(z["class"]) > > > On 07/24/2011 06:22 PM, Matevž Pavlič wrote: >> Hi list users, >> >> >> >> i am trying to do indicator kriging on a set of soild samples. >> >> >> >> I manage to do the predction of several classes into a >> "SpatialPixelsDataFrame" using the gstat object and predict() function. >> First a created a gstat object of several soil classes which I later fitted >> with a linear model of coregonalization. >> >> >> >> After that i used predict() to generate a SpatialPixelsDataFrame which hold >> information about predictions, variance and covariances. >> >> >> >> I am mostly interested in predoctions of each class. I know how to create a >> spplot or levelplot of each of the predicted classes on its own, but what i >> would like to achieve is to create a map of this predictions, so that the >> class with highest value (from 0-1) would prevail on each of the grid cell. >> In this way i would get a map of soil classes. >> >> >> >> Attached is the reproducible sample and a picture of prediction that i would >> like covert to a map of Soil classes (SoilFactor). But i get stuck when i >> try to create a actual map. That means getting the max value of each pixels >> and assign a name to that pixel. >> >> >> >> I hope i explained ok, if not i'll try again. >> >> >> >> Thanks for any ideas or help, m >> >> >> >> >> >> _______________________________________________ >> R-sig-Geo mailing list >> R-sig-Geo@r-project.org >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo > > -- > Edzer Pebesma > Institute for Geoinformatics (ifgi), University of Münster Weseler Straße > 253, 48151 Münster, Germany. Phone: +49 251 8333081, Fax: +49 251 8339763 > http://ifgi.uni-muenster.de > http://www.52north.org/geostatistics e.pebe...@wwu.de > > _______________________________________________ > R-sig-Geo mailing list > R-sig-Geo@r-project.org > https://stat.ethz.ch/mailman/listinfo/r-sig-geo -- Edzer Pebesma Institute for Geoinformatics (ifgi), University of Münster Weseler Straße 253, 48151 Münster, Germany. Phone: +49 251 8333081, Fax: +49 251 8339763 http://ifgi.uni-muenster.de http://www.52north.org/geostatistics e.pebe...@wwu.de _______________________________________________ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo