fdbk_dt_crps: Calculate CRPS(crps, crpsPot,Reli) from data.table applied on...

In rfxf/Rfdbk: Handling NetCDF feedback files

Description

Calculate CRPS(crps, crpsPot,Reli) from data.table applied on selected parts of the table (Caution, double check results! DT sorting might be modified!)

Usage

fdbk_dt_crps(DT, by)

Arguments

DT	data.table (columns 'veri_ens_member','obs' and 'veri_data' plus all variables to make forecasts distinguishable are required!!!)
by	stratify crps by (e.g. 'varno')

Value

data.table with columns as defined in 'by' plus scorename plus score

Author(s)

Felix <[e-mail]>

Examples

#EXAMPLE 1 (CRPS for each varno)
fnames                = system("ls ~/examplesRfdbk/talagrand/*SYNOP*",intern=T)
cond                  = list(veri_description="grepl('first guess ensemble member',veri_description)",
                             obs="!is.na(obs)",
                             statid="!is.na(as.numeric(statid)) & !duplicated(statid)",
                             veri_forecast_time="veri_forecast_time==100",
                             state="state%in%c(0,1,5)")
columnnames           = c("veri_data","varno","obs","veri_ens_member","veri_initial_date","veri_forecast_time","statid")
DT                    = fdbk_dt_multi_large(fnames,cond,columnnames,10)
DT[,varno:=varno_to_name(varno)]
fdbk_dt_crps(DT,by="varno")

#EXAMPLE 2 (CRPS decomosition for forecasts at SYNOP stations)
require(ggplot2)
fnames    = system("/bin/ls ~/examplesRfdbk/eps/*12/verSYNOP*",intern=T)
condition = list(veri_description="grepl('member',veri_description)",
                 state="state%in%c(0,1)",
                 statid="round(as.numeric(statid)/1000)==10 & !duplicated(statid)",
                 veri_forecast_time="veri_forecast_time>=1200")
columns   = c("obs","veri_data","varno","veri_ens_member","veri_forecast_time","statid","veri_initial_date")
DT        = fdbk_dt_multi_large(fnames,condition,columns,5)
CRPS = fdbk_dt_crps(DT,by=c("varno","veri_forecast_time"))
CRPS[,varno:=varno_to_name(varno,F)]
ggplot(CRPS,aes(x=veri_forecast_time,y=score))+geom_line()+geom_point()+facet_grid(~varno~scorename,scales="free_y")+theme_bw()+ggtitle("SYNOP (DE) CRPS decomosition")

#EXAMPLE 3 (slow...)(CRPS decomosition for european forecasts at TEMP stations)
require(ggplot2)
fnames    = system("/bin/ls ~/examplesRfdbk/eps/*12/verTEMP*",intern=T)
condition = list(veri_description="grepl('member',veri_description)",
                 state="state%in%c(0,1)",
                 level="level%in%c(100000,92500,85000,75000,70000,50000,40000,30000,25000,20000,10000)",
                 statid="round(as.numeric(statid)/1000)<=10 & !duplicated(statid)",
                 veri_forecast_time="veri_forecast_time>=1200",
                 varno="varno!=1")
columns   = c("obs","veri_data","varno","level","veri_ens_member","veri_forecast_time","veri_initial_date","statid")
DT        = fdbk_dt_multi_large(fnames,condition,columns,5)
CRPS      = fdbk_dt_crps(DT,by=c("varno","level","veri_forecast_time"))
CRPS[,varno:=varno_to_name(varno,F)]
ggplot(CRPS,aes(x=score,y=level,color=factor(veri_forecast_time),group=veri_forecast_time))+
 geom_path()+facet_wrap(~varno~scorename,scale="free_x",ncol=3)+
 scale_y_reverse()+theme_bw()+scale_colour_discrete("lead-time")

rfxf/Rfdbk documentation built on May 27, 2019, 7:22 a.m.