Re: [R-sig-Geo] Indicator kriging map creation
Hi Edzer, sorry for the long time to reply. >It might be annoying, and after all, it is your sole responsibility to make >sure the right labels end up in the right place. R can only prevent you from >making mistakes that are too obvious. What do you think would be the best thing to do about that? Predict the classes on 3D grid? I mean, i would like to create a code that would generate maps of several classes (variing numbers of them) each depth interval. I think that prediction on 3D grid would be the right thing? >Not so likely; I think it comes with simple kriging -- rare classes often show >weak autocorrelation, resulting in low estimates everywhere. >Without having seen your data, I suspect your classes 9 and 10 are like that. Attached are the co-variograms of ten classes. As i can see there are few classes (devet=level 10 , osem=level 9...) that have nugget (pure nugget effect). This is probably because there are only a few point with that class in the area : devet dva ena nic osem pet sedem set stiri tri 55085 209 3 783714989 Is there a way to fix this without discarding this classes? Regards, m <>___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] Indicator kriging map creation
On 07/24/2011 10:30 PM, Matevž Pavlič wrote:
> Hi,
>
> Thanks for the help! I realize I have a lot to learn , but sometimes you just
> get stuck.
>
> What i don't understad is this :
> when I use the first line of your code :
>> 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))
>
> and check levels :
>> levels(z$class)
> [1] "1" "2" "3" "4" "5" "6" "7" "8"
>
> the levels are numerical. I read in the help that " optional vector of the
> values that x might have taken. The default is the unique set of values taken
> by as.character(x), sorted into increasing order of x".
>
> I would like to create several IK maps with depth (every meter of depth).
> With doing that there will be areas that will not include all of the soil
> classes (or in this case factor levels). Will this be a problem with factor
> labels and map creating, if some of the soil classes will be missing?
It might be annoying, and after all, it is your sole responsibility to
make sure the right labels end up in the right place. R can only prevent
you from making mistakes that are too obvious.
>
> Also, you wrote " One or more of your classes never has the highest
> probability, it seems ".
> As i understand the apply functions it checks the rows in grid and assignes
> the maximum values of pred to the cell. So i do not understand how some
> predictions never have the highest probability since i have prediction in a
> SpatialGridDataFrame for each of the 10 classes? Could it be because the
> values in some cells are equal?
Not so likely; I think it comes with simple kriging -- rare classes
often show weak autocorrelation, resulting in low estimates everywhere.
Without having seen your data, I suspect your classes 9 and 10 are like
that.
Please do read ?which.max
>
> Thanks again for the answer and your time,
>
> matevz
>
>
> -----Original Message-
> From: Edzer Pebesma [mailto:edzer.pebe...@uni-muenster.de]
> Sent: Sunday, July 24, 2011 9:31 PM
> To: Matevž Pavlič
> Cc: r-sig-geo@r-project.org
> Subject: Re: [R-sig-Geo] Indicator kriging map creation
>
> 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
Re: [R-sig-Geo] Indicator kriging map creation
Hi,
Thanks for the help! I realize I have a lot to learn , but sometimes you just
get stuck.
What i don't understad is this :
when I use the first line of your code :
>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))
and check levels :
>levels(z$class)
[1] "1" "2" "3" "4" "5" "6" "7" "8"
the levels are numerical. I read in the help that " optional vector of the
values that x might have taken. The default is the unique set of values taken
by as.character(x), sorted into increasing order of x".
I would like to create several IK maps with depth (every meter of depth). With
doing that there will be areas that will not include all of the soil classes
(or in this case factor levels). Will this be a problem with factor labels and
map creating, if some of the soil classes will be missing?
Also, you wrote " One or more of your classes never has the highest
probability, it seems ".
As i understand the apply functions it checks the rows in grid and assignes the
maximum values of pred to the cell. So i do not understand how some
predictions never have the highest probability since i have prediction in a
SpatialGridDataFrame for each of the 10 classes? Could it be because the values
in some cells are equal?
Thanks again for the answer and your time,
matevz
-Original Message-
From: Edzer Pebesma [mailto:edzer.pebe...@uni-muenster.de]
Sent: Sunday, July 24, 2011 9:31 PM
To: Matevž Pavlič
Cc: r-sig-geo@r-project.org
Subject: Re: [R-sig-Geo] Indicator kriging map creation
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
Re: [R-sig-Geo] Indicator kriging map creation
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.
>>
>>
>>
&
Re: [R-sig-Geo] Indicator kriging map creation
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 25
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
