0

I wanted to computed zonal percentiles in R for raster values with a gridcode, using the DoBy package. Certain percentiles work perfectly, others produce more or less NA values, in this case -9999. What can be the cause of this discrepancy?

NOTE: All NA values appear in polygons that are not fully overlaping with the raster i'm interested in. This doesn't however explain why some percentile compute, others do not.

My applied percentile functions: 

quant1<-function(x,na.rm=TRUE){(quantile(x, c(.01), na.rm = TRUE))}
quant5<-function(x,na.rm=TRUE){(quantile(x, c(.05), na.rm = TRUE))}
quant15<-function(x,na.rm=TRUE){(quantile(x, c(.15), na.rm = TRUE))}
quant25<-function(x,na.rm=TRUE){(quantile(x, c(.25), na.rm = TRUE))}
quant75<-function(x,na.rm=TRUE){(quantile(x, c(.75), na.rm = TRUE))}
quant85<-function(x,na.rm=TRUE){(quantile(x, c(.85), na.rm = TRUE))}
quant95<-function(x,na.rm=TRUE){(quantile(x, c(.95), na.rm = TRUE))}
quant99<-function(x,na.rm=TRUE){(quantile(x, c(.99), na.rm = TRUE))}

Application of functions: 

q1<-summaryBy(RASTERVALU~grid_code,data=depth, FUN=quant1, keep.names=T)
q5<-summaryBy(RASTERVALU~grid_code,data=depth, FUN=quant5, keep.names=T)
q15<-summaryBy(RASTERVALU~grid_code,data=depth, FUN=quant15, keep.names=T)
q25<-summaryBy(RASTERVALU~grid_code,data=depth, FUN=quant25, keep.names=T)
q75<-summaryBy(RASTERVALU~grid_code,data=depth, FUN=quant75, keep.names=T)
q85<-summaryBy(RASTERVALU~grid_code,data=depth, FUN=quant85, keep.names=T)
q95<-summaryBy(RASTERVALU~grid_code,data=depth, FUN=quant95, keep.names=T)
q99<-summaryBy(RASTERVALU~grid_code,data=depth, FUN=quant99, keep.names=T)

This should to my understanding remove NA values before computation and produce results for all polygons, no matter the pixel count.

Follows an output table of the computation. Percentile 1 and 15 are most plagued by NA values.

q1  q5  q15 q25 q75 q85 q95 q99
-9999   -82.7934989929199   -9999   -9999   -82.6172466278076   -82.5455478668213   -82.4335483551025   -82.3479396057129
-9999   -82.0279998779297   -9999   -9999   -81.7740020751953   -81.6869964599609   -81.6319034576416   -81.6039762878418
-9999   -68.5245018005371   -9999   -68.6367511749268   -68.4617481231689   -68.4400024414062   -68.4128520965576   -68.318980255127
-9999   -68.1529998779297   -9999   -68.6242504119873   -67.9940013885498   -67.9098526000977   -67.8197998046875   -67.6305767059326
-9999   -66.0965003967285   -9999   -9999   -65.9724998474121   -65.9219970703125   -65.7605514526367   -65.6639614868164
-9999   -65.2804985046387   -9999   -9999   -64.8652515411377   -64.7608005523682   -64.5118511199951   -64.4869895172119
-9999   -62.7045001983643   -9999   -9999   -62.5639991760254   -62.4821016311646   -62.2751508712769   -62.1119895172119
-9999   -62.5089988708496   -9999   -62.6920013427734   -62.3634986877441   -62.2785495758057   -62.1639013290405   -62.1359615707397
-9999   -62.3390007019043   -9999   -9999   -62.1464996337891   -62.088650894165    -62.031251335144    -61.9182501983643
-9999   -58.6590003967285   -9999   -9999   -58.3725004196167   -58.1434007644653   -57.8779983520508   -57.6537393188477
-9999   -58.4934997558594   -9999   -58.8072500228882   -58.3682489395142   -58.3038509368896   -58.1767482757568   -58.0499892425537
-9999   -56.1360015869141   -9999   -56.3297491073608   -56.0120010375977   -55.9688512802124   -55.8899494171143   -55.8329501724243
-9999   -53.9230003356934   -9999   -54.2092504501343   -53.8269996643066   -53.7830009460449   -53.6978017807007   -53.6589603805542
-9999   -48.9314994812012   -9999   -9999   -48.7105007171631   -48.6123992919922   -48.5179513931274   -48.4739890289307
-9999   -42.2854995727539   -9999   -42.8615007400513   -41.8812494277954   -41.6038501739502   -41.3591510772705   -40.5169882965088
-9999   -42.1614990234375   -9999   -42.2462491989136   -42.0267505645752   -41.9099998474121   -41.8459506988525   -41.8219386291504
-9999   -38.1544990539551   -9999   -9999   -37.8874998092651   -37.6593996047974   -37.4445497512817   -37.3039208602905
-9999   -33.1689987182617   -9999   -33.2529983520508   -33.1027507781982   -33.0769996643066   -33.0347497940063   -33.0139999389648
-9999   -29.6655006408691   -9999   -9999   -29.5527496337891   -29.5060005187988   -29.4157503128052   -29.3299599266052
-9999   -29.638500213623    -9999   -29.7770004272461   -29.5520000457764   -29.5209999084473   -29.4968997955322   -29.4849806022644
-9999   -29.0570001602173   -9999   -9999   -28.9180006980896   -28.8618495941162   -28.7558505058289   -28.6939602470398
-9999   -27.1759996414185   -9999   -9999   -27.0627498626709   -27.0209999084473   -26.9890003204346   -26.9459791564941
-9999   -25.918999671936    -9999   -25.9869995117188   -25.8700008392334   -25.8515495300293   -25.8089500427246   -25.7189898681641

What could explain this effect?

Note: My data has no NA values. Certain grids just have a lower pixel count.

UPDATE: Yes, the data has -9999 for NA values, silly and easy mistake. Thanks for the help!

Castle
  • 53
  • 6

2 Answers2

1

From where I'm sitting, your result makes sense. -9999 are "classified" in the lower (0.1 and 0.15) quantiles.

> quantile(c(rep(-9999, 15), rnorm(50, mean = -50, sd = 5)), prob = c(0.1, 0.15, 0.25, 0.75))
        10%         15%         25%         75% 
-9999.00000 -9999.00000   -60.11869   -49.44009

Also note that according to your notation, quant5 function should probably be quant50.

Roman Luštrik
  • 69,533
  • 24
  • 154
  • 197
  • But what am i to make of this? I have values for all those pixels.. even if there's two pixles, i'd like to get a numeric result that the quantile 0.01 equals the lower value. And thank's for the heads up on the quant5-50 thing! – Castle Mar 23 '16 at 12:37
  • @Sander sorry, that wasn't part of your question (just a note). Please open a new question addressing this issue. Make sure you provide a [reproducible example](http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example). – Roman Luštrik Mar 23 '16 at 13:37
1

I would do:

library(raster)
r <- raster(ncols=10, nrows=10)
r[] <- runif(ncell(r)) * 1:ncell(r)
z <- r
z[] <- rep(1:5, each=20)

quants <- c(0.01, 0.05, 0.15, 0.25, 0.75, 0.85, 0.95, 0.99)

q <- zonal(r, z, function(i, ...) quantile(i, quants))
colnames(q) <- c("zone", quants)
q

Now, your data clearly has values of -9999. These are not NA, they are -9999. That is why your lowest quantile values have that value. If you want them to be NA, you may try:

NAvalue(r) <- -9999
q <- zonal(r, z, function(i, ...) quantile(i, quants, na.rm=TRUE))

or 

x <- reclassify(r, cbind(-9999, NA))
q <- zonal(r, z, function(i, ...) quantile(i, quants, na.rm=TRUE))
Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63