Jenks classification for raster representation

classic Classic list List threaded Threaded
3 messages Options
Reply | Threaded
Open this post in threaded view
|

Jenks classification for raster representation

Arnaud Mosnier
Dear list,

In order to allow others benefiting from my errors, see below a presentation
of a problem and it's solution by Roger Bivand (Thanks !).

######################################

I want to use natural breaks (jenks) method to find class intervals into
raster's values in order to plot it with an adapted color range.
I know that the package "classInt" allows to do exactly what I want but I do
not achieve to use it correctly.

Here is a small script that represent my case

# here a raster with some NAs

library(raster)
t<-raster(ncol = 10, nrow = 10)
t[1:2,] <- seq(1,10,1)
t[3:5,1:5] <- seq(1,5,1)
t[6:8,6:10] <-  seq(6,10,1)
t[9:10,] <- seq(1,10,1)
plot(t)

# Then I tried to use classIntervals and findColours function

library(classInt)
zClass <- classIntervals(as.data.frame(t)[!is.na(as.data.frame(t))], n=5,
style="jenks")  # I have to remove NAs
plot(t, col = findColours(zClass, pal = c("red", "yellow"))) # I do not
understand the result both in the plot and in the legend !

I would be pleased if you can give me tips about this.

########################################

I do not think that the plot method for RasterLayer objects lets you pass a
vector of breaks, nor does the image.plot from fields that maybe is the
source of the legend. If you do things simply, it works:

tSG <- as(t, "SpatialGridDataFrame")
spplot(tSG, "values", at=zClass$brks,
 col.regions=colorRampPalette(
c("red", "yellow"))(5))

or

image(tSG, "values", breaks=zClass$brks,
 col=colorRampPalette(c("red", "yellow"))(5))

or using graphics::image

image(as.image.SpatialGridDataFrame(tSG), breaks=zClass$brks,
 col=colorRampPalette(c("red", "yellow"))(5))

I'm a bit unsure about whether the itervals are closed which way.

It's up to the raster developers to answer.

Roger

--
Roger Bivand
Department of Economics, NHH Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: [hidden email]

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: Jenks classification for raster representation

Robert Hijmans
Arnaud,

You can (directly) do the same with Raster* objects what Roger did with sp
objects:

library(raster)
t<-raster(ncol = 10, nrow = 10)
t[1:2,] <- seq(1,10,1)
t[3:5,1:5] <- seq(1,5,1)
t[6:8,6:10] <-  seq(6,10,1)
t[9:10,] <- seq(1,10,1)
plot(t)

library(classInt)
zClass <- classIntervals(values(t), n=5,style="jenks")

# image
image(t, breaks=zClass$brks, col=colorRampPalette(c("red", "yellow"))(5) )

# or spplot
spplot(t, "values", at=zClass$brks, col.regions=colorRampPalette(c("red",
"yellow"))(5))

# or plot, after an adjustment to zClass that won't be needed in the future
(versions > 1.9-13)
zClass$brks[1] <- zClass$brks[1] - 1
plot(t, breaks=zClass$brks, col=colorRampPalette(c("red", "yellow"))(5) )


# The  only additional consideration is that classIntervals will fail, or be
very slow, for very large rasters.
# To avoid that, you can take a sample:

zClass <- classIntervals(na.omit(sampleRegular(t, 100000)),
n=5,style="jenks")
zClass$brks[1] <- zClass$brks[1] - 1

plot(t, breaks=zClass$brks, col=colorRampPalette(c("red", "yellow"))(5))


Robert


On Fri, Sep 16, 2011 at 12:58 PM, Arnaud Mosnier <[hidden email]>wrote:

> Dear list,
>
> In order to allow others benefiting from my errors, see below a
> presentation
> of a problem and it's solution by Roger Bivand (Thanks !).
>
> ######################################
>
> I want to use natural breaks (jenks) method to find class intervals into
> raster's values in order to plot it with an adapted color range.
> I know that the package "classInt" allows to do exactly what I want but I
> do
> not achieve to use it correctly.
>
> Here is a small script that represent my case
>
> # here a raster with some NAs
>
> library(raster)
> t<-raster(ncol = 10, nrow = 10)
> t[1:2,] <- seq(1,10,1)
> t[3:5,1:5] <- seq(1,5,1)
> t[6:8,6:10] <-  seq(6,10,1)
> t[9:10,] <- seq(1,10,1)
> plot(t)
>
> # Then I tried to use classIntervals and findColours function
>
> library(classInt)
> zClass <- classIntervals(as.data.frame(t)[!is.na(as.data.frame(t))], n=5,
> style="jenks")  # I have to remove NAs
> plot(t, col = findColours(zClass, pal = c("red", "yellow"))) # I do not
> understand the result both in the plot and in the legend !
>
> I would be pleased if you can give me tips about this.
>
> ########################################
>
> I do not think that the plot method for RasterLayer objects lets you pass a
> vector of breaks, nor does the image.plot from fields that maybe is the
> source of the legend. If you do things simply, it works:
>
> tSG <- as(t, "SpatialGridDataFrame")
> spplot(tSG, "values", at=zClass$brks,
>  col.regions=colorRampPalette(
> c("red", "yellow"))(5))
>
> or
>
> image(tSG, "values", breaks=zClass$brks,
>  col=colorRampPalette(c("red", "yellow"))(5))
>
> or using graphics::image
>
> image(as.image.SpatialGridDataFrame(tSG), breaks=zClass$brks,
>  col=colorRampPalette(c("red", "yellow"))(5))
>
> I'm a bit unsure about whether the itervals are closed which way.
>
> It's up to the raster developers to answer.
>
> Roger
>
> --
> Roger Bivand
> Department of Economics, NHH Norwegian School of Economics,
> Helleveien 30, N-5045 Bergen, Norway.
> voice: +47 55 95 93 55; fax +47 55 95 95 43
> e-mail: [hidden email]
>
>        [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: Jenks classification for raster representation

Tara Bridwell
Thank you Robert.  That worked perfectly.  I had been missing the
proper structure in using the rasterize function.

Tara

On Sat, Sep 17, 2011 at 5:42 PM, Robert J. Hijmans <[hidden email]> wrote:

> Arnaud,
>
> You can (directly) do the same with Raster* objects what Roger did with sp
> objects:
>
> library(raster)
> t<-raster(ncol = 10, nrow = 10)
> t[1:2,] <- seq(1,10,1)
> t[3:5,1:5] <- seq(1,5,1)
> t[6:8,6:10] <-  seq(6,10,1)
> t[9:10,] <- seq(1,10,1)
> plot(t)
>
> library(classInt)
> zClass <- classIntervals(values(t), n=5,style="jenks")
>
> # image
> image(t, breaks=zClass$brks, col=colorRampPalette(c("red", "yellow"))(5) )
>
> # or spplot
> spplot(t, "values", at=zClass$brks, col.regions=colorRampPalette(c("red",
> "yellow"))(5))
>
> # or plot, after an adjustment to zClass that won't be needed in the future
> (versions > 1.9-13)
> zClass$brks[1] <- zClass$brks[1] - 1
> plot(t, breaks=zClass$brks, col=colorRampPalette(c("red", "yellow"))(5) )
>
>
> # The  only additional consideration is that classIntervals will fail, or be
> very slow, for very large rasters.
> # To avoid that, you can take a sample:
>
> zClass <- classIntervals(na.omit(sampleRegular(t, 100000)),
> n=5,style="jenks")
> zClass$brks[1] <- zClass$brks[1] - 1
>
> plot(t, breaks=zClass$brks, col=colorRampPalette(c("red", "yellow"))(5))
>
>
> Robert
>
>
> On Fri, Sep 16, 2011 at 12:58 PM, Arnaud Mosnier <[hidden email]>wrote:
>
>> Dear list,
>>
>> In order to allow others benefiting from my errors, see below a
>> presentation
>> of a problem and it's solution by Roger Bivand (Thanks !).
>>
>> ######################################
>>
>> I want to use natural breaks (jenks) method to find class intervals into
>> raster's values in order to plot it with an adapted color range.
>> I know that the package "classInt" allows to do exactly what I want but I
>> do
>> not achieve to use it correctly.
>>
>> Here is a small script that represent my case
>>
>> # here a raster with some NAs
>>
>> library(raster)
>> t<-raster(ncol = 10, nrow = 10)
>> t[1:2,] <- seq(1,10,1)
>> t[3:5,1:5] <- seq(1,5,1)
>> t[6:8,6:10] <-  seq(6,10,1)
>> t[9:10,] <- seq(1,10,1)
>> plot(t)
>>
>> # Then I tried to use classIntervals and findColours function
>>
>> library(classInt)
>> zClass <- classIntervals(as.data.frame(t)[!is.na(as.data.frame(t))], n=5,
>> style="jenks")  # I have to remove NAs
>> plot(t, col = findColours(zClass, pal = c("red", "yellow"))) # I do not
>> understand the result both in the plot and in the legend !
>>
>> I would be pleased if you can give me tips about this.
>>
>> ########################################
>>
>> I do not think that the plot method for RasterLayer objects lets you pass a
>> vector of breaks, nor does the image.plot from fields that maybe is the
>> source of the legend. If you do things simply, it works:
>>
>> tSG <- as(t, "SpatialGridDataFrame")
>> spplot(tSG, "values", at=zClass$brks,
>>  col.regions=colorRampPalette(
>> c("red", "yellow"))(5))
>>
>> or
>>
>> image(tSG, "values", breaks=zClass$brks,
>>  col=colorRampPalette(c("red", "yellow"))(5))
>>
>> or using graphics::image
>>
>> image(as.image.SpatialGridDataFrame(tSG), breaks=zClass$brks,
>>  col=colorRampPalette(c("red", "yellow"))(5))
>>
>> I'm a bit unsure about whether the itervals are closed which way.
>>
>> It's up to the raster developers to answer.
>>
>> Roger
>>
>> --
>> Roger Bivand
>> Department of Economics, NHH Norwegian School of Economics,
>> Helleveien 30, N-5045 Bergen, Norway.
>> voice: +47 55 95 93 55; fax +47 55 95 95 43
>> e-mail: [hidden email]
>>
>>        [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
>        [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo