yuanchao-xu / hyfo Goto Github PK

Hello,

I have been trying to use the package “hyfo” for bias correction for my research project.

My research project aims to calculate drought projections for Scotland for the future using CMIP6 ensemble models. For this, I am trying to bias correct my CMIP6 climate model data using observational ERA5 datasets.

I have followed all the rules as outlined in the document but my code is not working and returns an error as follows :

Error in (function (classes, fdef, mtable) :
unable to find an inherited method for function ‘biasCorrect’ for signature ‘"RasterStackTS", "RasterStackTS", "RasterStackTS"’

I realise that there is something within my data that the function is not able to work with, hence the error, but I cannot figure out what it is.

I have attached a copy of my script (the part that covers the problem) and example data files that can be used to run the code.

If you can help me out, it will be immensely helpful !

Precipitation data attached :

1. ->historical climate model data

2. -> forecast climate model data

3. -> ERA5 observational dataset for the baseline period of study (1980-2010)

4. pr_Amon_BCC-CSM2-MR_historical_r1i1p1f1_gn_185001-201412_v20181126 2.nc

5. pr_Amon_BCC-CSM2-MR_ssp126_r1i1p1f1_gn_201501-210012_v20190314.nc

6. adaptor.mars.internal-1664806420.6580343-32648-12-01f0b7b0-a9fd-4fb3-b056-a49d9b5aa230 1.nc

biasFactor <- getBiasFactor(nc, nc1)
Error in if (!is.na(memberIndex)) { : argument is of length zero

Hello there,
I want to bais correct some data, my codes are following
new_data <- biasCorrect(ban_tas ,hindcast = ban_his, ban_obs, method="eqm")
But showing the following error

Error in biasCorrect.TS(frc, hindcast, obs, method, scaleType, preci, :
First column is not date or Wrong Date formate, check the format in ?as.Date{base}
and use as.Date to convert.If your input is a hyfo dataset, put input = "hyfo" as an
argument, check help for more info.

I have formated date to %Y-%m-%d and I also used as.Date function as instructed, but still showing the same error.
looking forward for your kind response.

How to convet my bias corrected output to netcdf file in your package by WRITEncdf ?
This is my output

newFrc11
$Variable
$Variable$varName
[1] "tas"

$xyCoords
$xyCoords$x
[1] 68.25 68.75 69.25 69.75 70.25 70.75 71.25 71.75 72.25 72.75 73.25 73.75 74.25 74.75 75.25 75.75 76.25 76.75 77.25 77.75
[21] 78.25 78.75 79.25 79.75 80.25 80.75 81.25 81.75 82.25 82.75 83.25 83.75 84.25 84.75 85.25 85.75 86.25 86.75 87.25 87.75
[41] 88.25 88.75 89.25 89.75 90.25 90.75 91.25 91.75 92.25 92.75 93.25 93.75 94.25 94.75 95.25 95.75 96.25 96.75 97.25
attr(,"name")
[1] "lon"

$xyCoords$y
[1] 7.75 8.25 8.75 9.25 9.75 10.25 10.75 11.25 11.75 12.25 12.75 13.25 13.75 14.25 14.75 15.25 15.75 16.25 16.75 17.25
[21] 17.75 18.25 18.75 19.25 19.75 20.25 20.75 21.25 21.75 22.25 22.75 23.25 23.75 24.25 24.75 25.25 25.75 26.25 26.75 27.25
[41] 27.75 28.25 28.75 29.25 29.75 30.25 30.75 31.25 31.75 32.25 32.75 33.25 33.75 34.25 34.75 35.25 35.75 36.25 36.75 37.25
attr(,"name")
[1] "lat"

$Dates
$Dates$start
[1] "1972-07-01" "1973-07-01" "1974-07-01" "1975-07-01" "1976-07-01" "1977-07-01" "1978-07-01" "1979-07-01" "1980-07-01"
[10] "1981-07-01" "1982-07-01" "1983-07-01" "1984-07-01" "1985-07-01" "1986-07-01" "1987-07-01" "1988-07-01" "1989-07-01"
[19] "1990-07-01" "1991-07-01" "1992-07-01" "1993-07-01" "1994-07-01" "1995-07-01" "1996-07-01" "1997-07-01" "1998-07-01"
[28] "1999-07-01" "2000-07-01"

$Data
, , 1

  [,1]     [,2]     [,3]     [,4]     [,5]     [,6]     [,7]     [,8]     [,9]    [,10]    [,11]    [,12]    [,13]

[1,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[2,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[3,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[4,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[5,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[6,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[7,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[8,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[9,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[10,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[11,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[12,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[13,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[14,] NA NA NA NA NA NA NA NA NA NA 25.29925 24.61353 24.67908
[15,] NA NA NA NA NA NA NA NA 25.12511 24.11722 23.53860 22.59216 23.05936
[16,] NA NA NA NA NA NA 25.77329 24.98765 23.75885 22.46367 22.24447 22.33635 23.07910
[,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24] [,25] [,26]
[1,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[2,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[3,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[4,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[5,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[6,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[7,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[8,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[9,] NA NA NA NA NA NA NA NA NA NA NA NA NA
[10,] NA NA NA NA NA NA NA NA 25.33255 24.88047 24.50799 24.39390 24.11259
[11,] NA NA NA NA 26.04569 25.45976 25.17771 24.91693 24.67589 24.60261 24.52710 24.46205 24.28529
[12,] NA NA 26.07548 25.52630 24.68124 24.11461 24.06714 23.95904 24.17682 24.83328 24.17924 24.71419 24.72445
[13,] 25.98223 25.32056 24.66288 24.14354 24.31315 24.10902 24.10460 24.20741 23.97003 24.14452 24.05985 24.22059 24.45213
[14,] 24.45948 23.81344 23.95213 24.35358 24.42159 24.47242 24.31926 24.30270 24.18323 24.43197 24.78346 24.27300 24.41734
[15,] 23.72778 24.16393 24.43164 24.68018 24.81718 24.72201 24.54731 24.55389 24.43554 24.41873 24.43338 24.47412 24.59300
[16,] 23.83839 24.42955 24.73748 24.94169 25.10924 25.17544 24.94215 24.86352 24.68728 24.49520 24.54827 24.67415 24.71652
[,27] [,28] [,29] [,30] [,31] [,32] [,33] [,34] [,35] [,36] [,37] [,38] [,39]
[1,] NA NA NA NA NA 24.13204 24.27227 NA NA NA NA NA NA
[2,] NA NA NA NA NA 24.17138 24.28241 24.47182 NA NA NA NA NA
[3,] NA NA 24.64331 24.50768 24.22869 26.23651 24.35168 24.51668 NA NA NA NA 23.99809
[4,] NA 24.94362 25.02834 25.03005 24.43764 24.54094 24.57757 24.51586 NA NA NA NA 24.10313
[5,] 25.29304 25.36167 25.31597 25.21395 25.05485 24.79874 24.84828 24.59642 NA NA 24.24702 24.06909 24.07320
[6,] 25.62296 25.68225 25.47631 25.36447 25.29600 25.04119 25.24482 24.92517 24.46755 24.23021 24.03012 24.09980 23.95408
[7,] 25.40863 25.72949 25.68777 25.59778 25.43232 25.34302 25.18136 25.03643 24.76089 24.17354 23.99676 23.93288 23.72056
[8,] 24.90822 25.70153 25.91480 25.83398 25.92430 25.32500 25.20195 25.01628 24.63886 24.27351 24.00937 23.58740 23.47781
[9,] NA 24.85315 25.38375 25.65846 26.22699 25.43312 25.20545 24.85020 24.48113 24.23078 24.00122 23.51918 23.41301
[10,] 23.72794 23.99383 25.19965 25.79520 26.20215 25.59667 24.92026 24.36112 24.07323 24.02584 23.95750 24.17184 23.34407
[11,] 24.24860 24.59470 25.07735 25.56223 25.63127 25.12458 24.66710 23.96674 23.63960 23.63918 23.69825 23.99610 23.43580
[12,] 24.39255 24.70620 24.98558 25.04645 25.13207 24.91902 24.31232 23.63337 23.44845 23.34163 23.42675 23.56369 23.54836
[13,] 24.52861 24.68137 25.04802 24.97884 24.73022 24.74731 24.47842 23.63171 23.25799 23.20684 23.22260 23.31351 23.30046
[14,] 25.09176 24.86698 24.95403 24.81645 24.48474 24.25193 23.80328 23.84610 23.46846 23.22048 23.26732 23.29955 23.29785
[15,] 24.91969 24.90235 25.01977 24.75814 24.48043 24.13700 23.55726 23.68934 23.57443 23.42813 23.43917 23.37698 23.31709
[16,] 24.87631 24.81226 24.99524 24.82615 24.67728 23.97280 23.54413 23.59028 23.32458 23.61204 23.61839 23.53930 23.47162
[,40] [,41] [,42] [,43] [,44] [,45] [,46] [,47] [,48] [,49] [,50]
[1,] NA NA NA NA NA NA NA NA NA NA NA
[2,] NA NA NA NA NA NA NA NA NA NA NA
[3,] NA NA NA NA NA NA NA NA NA NA NA
[4,] 27.904139 28.057581 NA NA NA NA NA NA NA NA NA
[5,] 24.066264 23.856365 NA NA NA NA NA NA NA NA NA
[6,] 23.945708 23.792701 NA NA NA NA NA NA NA NA NA
[7,] 23.695852 23.557860 NA NA NA NA NA NA NA NA NA
[8,] 23.501514 23.510431 23.687422 NA NA NA NA NA NA NA NA
[9,] 23.336036 23.253902 23.548784 23.647974 NA NA NA NA NA NA NA
[10,] 23.129014 23.148600 23.388872 23.532022 23.5719471 NA NA NA NA NA NA
[11,] 23.115440 23.226533 23.246485 23.360163 23.4999485 23.531656 NA NA NA NA NA
[12,] 23.117430 23.109285 23.102739 23.139566 23.4932434 23.772039 23.49490656 NA NA NA NA
[13,] 23.199348 23.081411 23.330530 23.148776 23.2490951 23.410297 23.29084315 22.882667 NA NA NA
[14,] 23.237298 23.160778 23.238129 23.009823 23.1568906 23.356963 23.08178030 22.553595 21.935429 20.6813363 18.457172
[15,] 23.306511 23.196395 23.045060 22.963134 23.1215545 23.089809 22.64315224 22.122121 21.183299 19.4849491 16.443563
[16,] 23.269439 23.240685 23.141365 23.062956 23.1732535 22.787622 22.26169995 21.538460 19.854678 17.1779172 13.479879
[,51] [,52] [,53] [,54] [,55] [,56] [,57] [,58] [,59] [,60]
[1,] NA NA NA NA NA NA NA NA NA NA
[2,] NA NA NA NA NA NA NA NA NA NA
[3,] NA NA NA NA NA NA NA NA NA NA
[4,] NA NA NA NA NA NA NA NA NA NA
[5,] NA NA NA NA NA NA NA NA NA NA
[6,] NA NA NA NA NA NA NA NA NA NA
[7,] NA NA NA NA NA NA NA NA NA NA
[8,] NA NA NA NA NA NA NA NA NA NA
[9,] NA NA NA NA NA NA NA NA NA NA
[10,] NA NA NA NA NA NA 0 0 NA NA
[11,] NA NA NA 10.065542 NA NA 0 0 0 NA
[12,] NA 15.859998 11.853418 7.751617 3.1953096 0 0 0 0 NA
[13,] 17.236752 13.797530 9.487059 5.098226 0.9220651 0 0 0 0 NA
[14,] 15.349826 11.037953 6.607987 2.261481 0.0000000 0 0 0 0 NA
[15,] 12.473914 7.832510 3.209839 0.000000 0.0000000 0 0 0 0 NA
[16,] 8.540425 3.694814 0.000000 0.000000 0.0000000 0 0 0 0 NA

[ reached getOption("max.print") -- omitted 43 row(s) and 28 matrix slice(s) ]
attr(,"dimensions")
[1] "lon" "lat" "time"

$Loaded
[1] "by hyfo package, https://yuanchao-xu.github.io/hyfo/"

$biasCorrected_by
[1] "eqm"

Hello,
I am trying to perform the {R-Hyfo} bias correction of total monthly precipitation data from the ETA-HadGEM2-ES, ETA-Miroc5 and ENSEMBLE regionalized climate models, based on specialized observed data called “Xavier” (DOI: 10.1002/joc.4518).

However, when I run the code below, I get the following error:
Error in biasCorrect_core(frcData[lon, lat, ], hindcastData[lon, lat, : This method needs frc data to have the same length as obs data.

The same data used in the correction with the Delta method did not obtain errors when used for the correction with the EQM, GQM, SCALING methods; and I can't identify what is wrong.

My question is: Can this method only correct the bias if they have the same number of data for observations and frc? Also, is the default method 'delta' (Fig. 1) or 'scaling' (Fig. 2)?

Fig 1 https://github.com/Yuanchao-Xu/manuals/blob/master/hyfo.Rmd (Line 783)
Fig 2 https://cran.r-project.org/web/packages/hyfo/hyfo.pdf (Page 7)

I am sending below my code and data. Could you please check it?
Data.zip

Figure 1

Figure 2

Thanks.

My code:

#------------------------------------ Script -------------------------------------------------------
V.select.month <- c("JAN", "FEB", "MAR", "APR", "MAY", "JUN", "JUL", "AUG", "SEP", "OCT", "NOV", "DEC")

for(i in 1:12) { # Main loop is initiated

#---------------------- Selection for each month
#1980 - 2000
Prec_hist_1980_2000_Xavier <- hyfo::loadNcdf(filePath="./Prec_monthly_historical_Xavier_1980-2000.nc", varname="Precipitation", year = 1980:2000, month = i)
Prec_hist_1980_2000_HadGEM2 <- hyfo::loadNcdf(filePath="./Prec_monthly_historical_ETA-HadGEM2-ES_1980-2000_Doce.nc", varname="Precipitation", year = 1980:2000, month = i)
Prec_hist_1980_2000_MIROC <- hyfo::loadNcdf(filePath="./Prec_monthly_historical_ETA-Miroc5_1980-2000_Doce.nc", varname="Precipitation", year = 1980:2000, month = i)
Prec_hist_1980_2000_ENSEMBLE <- hyfo::loadNcdf(filePath="./Prec_monthly_historical_ENSEMBLE_1980-2000_Doce.nc", varname="Precipitation", year = 1980:2000, month = i)

#2001 - 2005
#Prec_hist_2001_2005_Xavier <- hyfo::loadNcdf(filePath= "./Prec_monthly_historical_Xavier_2001-2005.nc", varname="Precipitation", year = 2001:2005, month = i)
Prec_hist_2001_2005_HadGEM2 <- hyfo::loadNcdf(filePath="./Prec_monthly_historical_ETA-HadGEM2-ES_2001-2005_Doce.nc", varname="Precipitation", year = 2001:2005, month = i)
Prec_hist_2001_2005_MIROC <- hyfo::loadNcdf(filePath="./Prec_monthly_historical_ETA-Miroc5_2001-2005_Doce.nc", varname="Precipitation", year = 2001:2005, month = i)
Prec_hist_2001_2005_ENSEMBLE <- hyfo::loadNcdf(filePath= "./Prec_monthly_historical_ENSEMBLE_2001-2005_Doce.nc", varname="Precipitation", year = 2001:2005, month = i)

#---------------------- Correction
# Delta biasCorrect {hyfo} is used to bias-correct the input timeseries hyfo dataset for the period 2001 a 2005
HadGEM2_2001_2005.DELTA.BC <- hyfo::biasCorrect(obs = Prec_hist_1980_2000_Xavier, hindcast = Prec_hist_1980_2000_HadGEM2, frc = Prec_hist_2001_2005_HadGEM2, method = 'delta', preci = TRUE)
MIROC_2001_2005.DELTA.BC <- hyfo::biasCorrect(obs = Prec_hist_1980_2000_Xavier, hindcast = Prec_hist_1980_2000_MIROC, frc = Prec_hist_2001_2005_MIROC, method = 'delta', preci = TRUE)
ENSEMBLE_2001_2005.DELTA.BC <- hyfo::biasCorrect(obs = Prec_hist_1980_2000_Xavier, hindcast = Prec_hist_1980_2000_ENSEMBLE, frc = Prec_hist_2001_2005_ENSEMBLE, method = 'delta', preci = TRUE)

#---------------------- Export
#DELTA writeNcdf {hyfo} is used write raster objects to netCDF files for the period 2001 a 2005
writeNcdf(gridData=HadGEM2_2001_2005.DELTA.BC,filePath=paste0("HadGEM2_2001_2005_DELTA_BC_",V.select.month[i],".nc"))
writeNcdf(gridData=MIROC_2001_2005.DELTA.BC,filePath=paste0("MIROC_2001_2005_DELTA_BC_",V.select.month[i],".nc"))
writeNcdf(gridData=ENSEMBLE_2001_2005.DELTA.BC,filePath=paste0("ENSEMBLE_2001_2005_DELTA_BC_",V.select.month[i],".nc"))

} #Main loop is closed

#END OF CODE

You will be aware, for example from:
https://r-spatial.org/r/2022/04/12/evolution.html,
https://r-spatial.org/r/2022/12/14/evolution2.html,
https://r-spatial.org/r/2023/04/10/evolution3.html and
https://rsbivand.github.io/csds_jan23/bivand_csds_ssg_230117.pdf and
perhaps view https://www.youtube.com/watch?v=TlpjIqTPMCA&list=PLzREt6r1NenmWEidssmLm-VO_YmAh4pq9&index=1
that rgdal, rgeos and maptools will be retired this
year, in October 2023.

'rgdal', 'maptools', 'rgeos' are declared as strong dependencies, but only rgdal::readOGR is found by pkgapi, NAMESPACE imports from all three, but rgdal::readOGR is only used in R/shp2cat.R, line 49. See https://r-spatial.org/r/2023/04/10/evolution3.html for alternatives terra::vect or sf::st_read, and any other missing functions. Run CMD check in a scenario using sp evolution status 2 (substitute use of rgdal with sf for projection/transformation/CRS) and absence of retiring packages from the library path. Fix best by June, latest by October 2023.

Dear Yuanchao,

I have a daily time series of precipitation, which I intend to transform it into a single NetCDF file of daily precipitation. So, first, I would like to explore the best way to impute missing values, grid the daily values of the different available weather stations and get daily raster files into a .nc format. I have learned about the Package "hyfo" and I think I can do it through it. However, I am experiencing some issues, and thus I haven't found a straightforward way to do it.

My database contains id, date, cod, lon, lat, elev, and precipitation value of 41 climate stations for the period 1980-2018. I have NAs values in the database. I organised the dataset following the hyfo Easy Start, so the first column presents date and the other 50 columns present a climate station, while rows present the values/NAs [14245 observations of 42 variables], as shown in the image below. The location of the stations [long, lat, elev, cod] is kept in a different file.

The annual precipitation values and NA values per each station look like this:

Since I have missing values, I attempted to use the 'FillGap' function, but I got the following error message: "Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...) : 0 (non-NA) cases" what I understand is that some stations with missing values are making noise so the function gets the error. Not so sure whether or not deleting those station would allow me to fill gaps with that function because I already deleted others. I am also not so sure that a method based on the correlation between neighbouring stations is the best approach to fill data gaps as the stations are located at a very irregular relieve. So my first question is what might be wrong when filling gaps? There might be a better function?

On the other hand, I learned that Package hyfo presents functions to read, analyse and write raster and NetCDF forms. However, I do not find a way out to grid the daily time series and get the .nc file. I would appreciate any input, idea or suggestion.

Thank you,

Sincerely,

Sir,

Please help to fix this problem.

filePath <- system.file("extdata", "uas20062050rcp45concat-cdo.nc", package = "hyfo")
uas <- getNcdfVar(filePath)
nc <- loadNcdf(filePath, uas)
Error in loadNcdf(filePath, uas) : Your file has less than 3 dimensions.

regards,
Subhasis

Hi,

I am trying to bias correct for time series data and I get this error:
Error in biasCorrect.TS(frc, hindcast, obs, method, scaleType, preci, :
First column is not date or Wrong Date formate, check the format in ?as.Date{base}
and use as.Date to convert.If your input is a hyfo dataset, put input = "hyfo" as an
argument, check help for more info.

This is the line of code I am using to bias correct:
PREQM26 <- biasCorrect(forecastP26 ,hindcast = hindcastP26, TRMM_DF,method="eqm")

where forecastP26 is my forecast data, hindcastP26 is my hindcast data and TRMM_DF is my observed data. I double check that each one of these dataframes does have the first column as the date and use the as.Date function to convert the first column to the correct date per the instructions of the error. After all of that, I still have received the same error. Any and all feedback on this issue is greatly appreciated!

Hi
I have attached a sample file in which I am getting the following error:

filePath <- system.file("extdata", "tasmax_WAS-44_ECMWF-ERAINT_evaluation_r1i1p1_IITM-RegCM4-4_v5_day_19790101-19801231.nc", package = "hyfo")
varname <- getNcdfVar(filePath)
ev <- loadNcdf(filePath, varname)
Error in loadNcdf(filePath, varname) :
Your file has less than 3 dimensions.

As my data contain three dimensions in it.
Wait for a reply.

I recently used you Hyfo package in order to use the bias correction part for my time series data. I faced with data format error(Error in as.POSIXlt.character(as.character(x), ...) :
character string is not in a standard unambiguous format), though i made my time series data in excel first column date and second column the value. would you please help me a little how to solve this problem.
additional information
I used bias correction function and my time series data are monthly available.
Thanks
Elnaz,

I have one dataset with 33 stations for observed station data and hindcast for the same stations, so I need to get the bias factor for each station and apply it to correct the projection for each station.I have been following the guideline but when i input my dataset,i get an error;" Error in getClass(Class, where = topenv(parent.frame())) :
“biasFactor.multiMember” is not a defined class".

T2_raj_2013-01-04.zip
Hi
When I am running the following command I am again getting the same error, irrespective of I changed the varname=Tas i,e name of the variable.

filePath <- system.file("extdata", "Tas_raj_2013-01-04.nc", package = "hyfo")
varname <- getNcdfVar(filePath)
histrocial <- loadNcdf(filePath, "Tas")
Error in loadNcdf(filePath, "Tas") : Your file has less than 3 dimensions.

Please help me solve this problem. I am attaching my sample file below.
Tas_raj_2013-01-04.zip