That's true, but I get as result one value which is the same for all the raster points.
I thought I could apply this function to each point by using the [] but it doesn't work...
for example, if after 

library(verification)

library(raster)

set.seed(1)

r <- raster(ncol=10, nrow=10)

pred1 <- init(r, fun=runif)

pred2 <- init(r, fun=runif)

pred3 <- init(r, fun=runif)

EUE <- stack(pred1, pred2, pred3)

WB_LT <- init(r, fun=runif)
EUE[1]<-NA

WB_LT[1]<-NA


# you run:



CRPS2 <- raster(WB_LT)
CRPS2[]<-crpsDecomposition(na.omit(WB_LT[]),na.omit(EUE[]))$CRPS
CRPS2[2]

# you get a different result from 

crpsDecomposition(na.omit(WB_LT[2]),na.omit(EUE[2]))$CRPS  

#The correct result could be obtained with
for (i in 1:100){
   CRPS2[i]<-crpsDecomposition(na.omit(WB_LT[i]),na.omit(EUE[i]))$CRPS  
}

#but I'd like to avoid the for loop as not efficient
Also, the na.omit doesn't overcome the problem in this case, so it crashes for points containing NA. In the example above it only runs if you start from i=2 (as i=1 has NA)

for (i in 2:100){

   CRPS2[i]<-crpsDecomposition(na.omit(WB_LT[i]),na.omit(EUE[i]))$CRPS  

}


Any suggestion?

Lorenzo




Date: Thu, 5 Sep 2013 08:50:15 -0700
Subject: Re: [R-sig-Geo] calculate CRPS on rasters
From: vitale232@gmail.com
To: alfios17@hotmail.com
CC: r-sig-geo@r-project.org

This seems to run:CRPS2[]<-crpsDecomposition(na.omit(WB_LT[]),na.omit(EUE[]))$CRPS


On Thu, Sep 5, 2013 at 6:43 AM, Lorenzo Alfieri <alfios17@hotmail.com> wrote:

Dear list,

following this thread, I'm having another problem.

I'm now using the command crpsDecomposition, of the verification package.

It's quite similar to the crps.





library(verification)

library(raster)

set.seed(1)

r <- raster(ncol=10, nrow=10)

pred1 <- init(r, fun=runif)

pred2 <- init(r, fun=runif)

pred3 <- init(r, fun=runif)

EUE <- stack(pred1, pred2, pred3)

WB_LT <- init(r, fun=runif)

CRPS <- raster(WB_LT)

CRPS[]<-crpsDecomposition(WB_LT[],EUE[])$CRPS



Up to here it works fine.

Problems start when there are NA values  in the input.

For example:



EUE[1]<-NA

WB_LT[1]<-NA

CRPS[]<-crpsDecomposition(WB_LT[],EUE[])$CRPS



In this case I get as output:



Error in prev[index, 1] - obs[index] :

  non-numeric argument to binary operator



In fact if I run:

crpsDecomposition(WB_LT[1],EUE[1])$CRPS

I get as output:

  Error in which(obs < prev[, 1]) : subscript out of bounds



One option is replacing all NAs with zeros

EUE[is.na(EUE)]<-0

before running

CRPS[]<-crpsDecomposition(WB_LT[],EUE[])$CRPS

but it's not very correct and in addition it takes ages, as the EUE is a big variable in my case.



Anyone can help me address this issue?



Regards



Lorenzo









> From: alfios17@hotmail.com

> To: r.hijmans@gmail.com

> CC: r-sig-geo@r-project.org

> Subject: RE: [R-sig-Geo] calculate CRPS on rasters

> Date: Fri, 30 Nov 2012 09:49:00 +0100

>

>

> Robert,

> you're right, it works.

> I must have done something wrong yesterday

> I'm trying to optimize the computation time

> I've tested the two alternatives on a predictor of 51 layers of 680x810 pixels

>

> CRPS[] <- crps(obs[],pred[])$crps                 took 0.98 second

> CRPS  <- overlay(obs, pred, fun=function(x,y) crps(x,y)$crps)          took 1.45 seconds

>

> anyway, thank you all!

>

> Lorenzo

>

>

> ________________________________

> > Date: Thu, 29 Nov 2012 10:34:36 -0800

> > Subject: Re: [R-sig-Geo] calculate CRPS on rasters

> > From: r.hijmans@gmail.com

> > To: alfios17@hotmail.com

> > CC: r-sig-geo@r-project.org

> >

> > Lorenzo,

> >

> > The below (your example) works for me, using the overlay as I suggested.

> > Perhaps you did not load the verification package?

> >

> > Robert

> >

> > library(verification)      #includes function crps

> > library(raster)

> > set.seed(1)

> > r <- raster(ncol=10, nrow=10)

> > pred1 <- init(r, fun=runif)

> > pred2 <- init(r, fun=runif)

> > pred3 <- init(r, fun=runif)

> > pred <- stack(pred1, pred2, pred3)

> > obs <- init(r, fun=runif)

> >

> > CRPS  <- overlay(obs, pred, fun=function(x,y) crps(x,y)$crps)

> >

> >

> >  > CRPS

> > class       : RasterLayer

> > dimensions  : 10, 10, 100  (nrow, ncol, ncell)

> > resolution  : 36, 18  (x, y)

> > extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)

> > coord. ref. : +proj=longlat +datum=WGS84

> > data source : in memory

> > names       : layer

> > values      : 0.003277372, 0.7726679  (min, max)

> >

> >

> > # More fancy, to get all the crps output:

> >

> > x  <- overlay(obs, pred, fun=function(x,y) as.matrix(data.frame(crps(x,y))))

> >

> >

> >

> > On Thu, Nov 29, 2012 at 1:15 AM, Lorenzo Alfieri

> > <alfios17@hotmail.com<mailto:alfios17@hotmail.com>> wrote:

> >

> > Hi Robert,

> > I was hoping to use the overlay function too, but it seems it doesn't

> > work in this case. See the result:

> >

> > CRPS  <- overlay(obs, pred, fun=function(x,y) crps(x,y)$crps)

> >

> > Error in function (x, fun, filename = "", recycle = TRUE, ...)  :

> >    cannot use this formula, probably because it is not vectorized

> >

> > Lorenzo

> >

> >

> >

> > Date: Wed, 28 Nov 2012 11:25:36 -0800

> > Subject: Re: [R-sig-Geo] calculate CRPS on rasters

> > From: r.hijmans@gmail.com<mailto:r.hijmans@gmail.com>

> > To: alfios17@hotmail.com<mailto:alfios17@hotmail.com>

> > CC: etiennebr@gmail.com<mailto:etiennebr@gmail.com>;

> > r-sig-geo@r-project.org<mailto:r-sig-geo@r-project.org>

> >

> > Lorenzo,

> > Given that solution, I think you can also express this (in a

> > memory-safe fashion) like this :

> > CRPS  <- overlay(obs, pred, fun=function(x,y) crps(x, y)$crps)

> >

> > Robert

> >

> > On Tue, Nov 27, 2012 at 3:11 AM, Lorenzo Alfieri

> > <alfios17@hotmail.com<mailto:alfios17@hotmail.com>> wrote:

> >

> >

> >

> >

> > Etienne,

> >

> > thank you for the tip

> >

> > Now it runs trough, by using

> >

> >

> >

> > CRPS <- raster(ncol=10, nrow=10)

> >

> > CRPS[] <- crps(obs[],pred[])$crps

> >

> >

> >

> > Lorenzo

> >

> >

> >

> >

> >

> > Date: Mon, 26 Nov 2012 21:37:13 -0500

> >

> > Subject: Re: [R-sig-Geo] calculate CRPS on rasters

> >

> > From: etiennebr@gmail.com<mailto:etiennebr@gmail.com>

> >

> > To: alfios17@hotmail.com<mailto:alfios17@hotmail.com>

> >

> > CC: r-sig-geo@r-project.org<mailto:r-sig-geo@r-project.org>

> >

> >

> >

> > Lorenzo,

> >

> >

> >

> > I don't know about your specific function, but you can access matrix

> > using the []'s. So maybe something like :

> >

> > CRPS <- raster(ncol=10, nrow=10)

> >

> > CRPS[] <- crps(obs[],pred[])

> >

> >

> >

> >

> >

> > Etienne

> >

> >

> >

> > 2012/11/26 Lorenzo Alfieri

> > <alfios17@hotmail.com<mailto:alfios17@hotmail.com>>

> >

> >

> >

> > CRPS <- raster(ncol=10, nrow=10)

> >

> >

> >

> > for (i in 1:length(obs)){

> >

> >

> >

> >    CRPS[i] <- crps(obs[i],pred[i])$CRPS

> >

> >

> >

> >          [[alternative HTML version deleted]]

> >

> >

> >

> > _______________________________________________

> >

> > R-sig-Geo mailing list

> >

> > R-sig-Geo@r-project.org<mailto:R-sig-Geo@r-project.org>

> >

> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo

> >

> >

> >

> >

>



        [[alternative HTML version deleted]]



_______________________________________________

R-sig-Geo mailing list

R-sig-Geo@r-project.org

https://stat.ethz.ch/mailman/listinfo/r-sig-geo



-- 
Andrew P. VitaleMasters StudentDepartment of GeographyUniversity of Nevada, Reno(412) 915-3632
vitale232@gmail.com
 		 	   		  
	[[alternative HTML version deleted]]

