Skip to content

Commit a102554

Browse files
added multiplicative procedure for QM and QDM
1 parent fecdd7f commit a102554

File tree

2 files changed

+63
-34
lines changed

2 files changed

+63
-34
lines changed

cmethods/CMethods.py

Lines changed: 40 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -423,7 +423,8 @@ def quantile_mapping(cls,
423423
Add (+):
424424
(1.) X^{*QM}_{sim,p}(i) = F^{-1}_{obs,h} \left\{F_{sim,h}\left[X_{sim,p}(i)\right]\right\}
425425
Mult (*):
426-
maybe the same
426+
(1.) --//--
427+
427428
----- R E F E R E N C E S -----
428429
429430
Multiplicative implementeation by 'deleted profile' OR Adrian Tompkins [email protected], posted on November 8, 2016 at
@@ -442,15 +443,12 @@ def quantile_mapping(cls,
442443

443444
cdf_obs = cls.get_cdf(obs, xbins)
444445
cdf_simh = cls.get_cdf(simh, xbins)
445-
epsilon = np.interp(simp, xbins, cdf_simh) # Eq. 1
446-
446+
epsilon = np.interp(simp, xbins, cdf_simh) # Eq. 1
447447
return cls.get_inverse_of_cdf(cdf_obs, epsilon, xbins) # Eq. 1
448-
448+
449449
elif kind == '*':
450450
''' Inspired by Adrian Tompkins [email protected] posted here:
451451
https://www.researchgate.net/post/Does-anyone-know-about-bias-correction-and-quantile-mapping-in-PYTHON
452-
- note that values exceeding the range of the training set
453-
are set to -999 at the moment - possibly could leave unchanged?
454452
'''
455453

456454
obs, simh, simp = np.sort(obs), np.sort(simh), np.array(simp)
@@ -545,27 +543,38 @@ def quantile_delta_mapping(cls,
545543
------ E Q U A T I O N S -----
546544
547545
Add (+):
548-
(1) \varepsilon(i) = F_{sim,p}\left[X_{sim,p}(i)\right], \hspace{1em} \varepsilon(i)\in\{0,1\}
546+
(1.1) \varepsilon(i) = F_{sim,p}\left[X_{sim,p}(i)\right], \hspace{1em} \varepsilon(i)\in\{0,1\}
549547
550-
(2) X^{QDM(1)}_{sim,p}(i) = F^{-1}_{obs,h}\left[\varepsilon(i)\right]
548+
(1.2) X^{QDM(1)}_{sim,p}(i) = F^{-1}_{obs,h}\left[\varepsilon(i)\right]
551549
552-
(3) \Delta(i) & = F^{-1}_{sim,p}\left[\varepsilon(i)\right] - F^{-1}_{sim,h}\left[\varepsilon(i)\right] \\[1pt]
550+
(1.3) \Delta(i) & = F^{-1}_{sim,p}\left[\varepsilon(i)\right] - F^{-1}_{sim,h}\left[\varepsilon(i)\right] \\[1pt]
553551
& = X_{sim,p}(i) - F^{-1}_{sim,h}\left\{F^{}_{sim,p}\left[X_{sim,p}(i)\right]\right\}
554552
555-
(4) X^{*QDM}_{sim,p}(i) = X^{QDM(1)}_{sim,p}(i) + \Delta(i)
553+
(1.4) X^{*QDM}_{sim,p}(i) = X^{QDM(1)}_{sim,p}(i) + \Delta(i)
556554
557555
Mult (*):
556+
(1.1) --//--
557+
558+
(1.2) --//--
559+
560+
(2.3) \Delta(i) & = \frac{ F^{-1}_{sim,p}\left[\varepsilon(i)\right] }{ F^{-1}_{sim,h}\left[\varepsilon(i)\right] } \\[1pt]
561+
& = \frac{ X_{sim,p}(i) }{ F^{-1}_{sim,h}\left\{F^{}_{sim,p}\left[X_{sim,p}(i)\right]\right\} }
562+
563+
(2.4) X^{*QDM}_{sim,p}(i) = X^{QDM(1)}_{sim,p}(i) \cdot \Delta(i)
558564
559565
----- R E F E R E N C E S -----
560566
Tong, Y., Gao, X., Han, Z. et al. Bias correction of temperature and precipitation over China for RCM simulations using the QM and QDM methods. Clim Dyn 57, 1425–1443 (2021).
561567
https://doi.org/10.1007/s00382-020-05447-4
568+
569+
----- N O T E S -----
570+
571+
@param global_min: float | this parameter can be set when kind == '*' to define a custom lower limit. Otherwise 0.0 is used.
562572
563573
'''
564574

565575
if group != None: return cls.grouped_correction(method='quantile_delta_mapping', obs=obs, simh=simh, simp=simp, group=group, n_quantiles=n_quantiles, kind=kind, **kwargs)
566576
elif kind == '+':
567577
obs, simh, simp = np.array(obs), np.array(simh), np.array(simp) # to achieve higher accuracy
568-
569578
global_max = max(np.amax(obs), np.amax(simh))
570579
global_min = min(np.amin(obs), np.amin(simh))
571580
wide = abs(global_max - global_min) / n_quantiles
@@ -576,14 +585,27 @@ def quantile_delta_mapping(cls,
576585
cdf_simp = cls.get_cdf(simp, xbins)
577586

578587
# calculate exact cdf values of $F_{sim,p}[T_{sim,p}(t)]$
579-
epsilon = np.interp(simp, xbins, cdf_simp) # Eq. 1
580-
581-
QDM1 = cls.get_inverse_of_cdf(cdf_obs, epsilon, xbins) # Eq. 2
582-
delta = simp - cls.get_inverse_of_cdf(cdf_simh, epsilon, xbins) # Eq. 3
588+
epsilon = np.interp(simp, xbins, cdf_simp) # Eq. 1.1
589+
QDM1 = cls.get_inverse_of_cdf(cdf_obs, epsilon, xbins) # Eq. 1.2
590+
delta = simp - cls.get_inverse_of_cdf(cdf_simh, epsilon, xbins) # Eq. 1.3
591+
return QDM1 + delta # Eq. 1.4
592+
593+
elif kind == '*':
594+
obs, simh, simp = np.array(obs), np.array(simh), np.array(simp)
595+
global_max = max(np.amax(obs), np.amax(simh))
596+
wide = global_max / n_quantiles
597+
xbins = np.arange(kwargs.get('global_min', 0.0), global_max + wide, wide)
583598

584-
return QDM1 + delta # Eq. 4
585-
elif kind == '*': raise ValueError('Only "+" is implemented!')
586-
else: raise ValueError('Only "+" is implemented!')
599+
cdf_obs = cls.get_cdf(obs,xbins)
600+
cdf_simh = cls.get_cdf(simh,xbins)
601+
cdf_simp = cls.get_cdf(simp, xbins)
602+
603+
epsilon = np.interp(simp, xbins, cdf_simp) # Eq. 1.1
604+
QDM1 = cls.get_inverse_of_cdf(cdf_obs, epsilon, xbins) # Eq. 1.2
605+
delta = simp / cls.get_inverse_of_cdf(cdf_simh, epsilon, xbins) # Eq. 2.3
606+
return QDM1 * delta # Eq. 2.4
607+
608+
else: raise ValueError(f'Unknown kind {kind}!')
587609

588610
# * -----========= G E N E R A L =========------
589611
@staticmethod

0 commit comments

Comments
 (0)