Statistical Post Processing

1. Statistical Post Processing Hong Guan Bo Cui and Yuejian Zhu EMC/NCEP/NWS/NOAA Presents for NWP Forecast Training Class March 31, 2015, Fuzhou, Fujian, China

2. 1. Background

3. North American Ensemble Forecast System (NAEFS) International project to produce operational multi-center ensemble products Bias correction and combines global ensemble forecasts from Canada & USA Generates products for: Weather forecasters Specialized users End users Operational outlet for THORPEX research using TIGGE archive

4. Statement The North American Ensemble Forecast System (NAEFS) combines state of the art weather forecast tools, called ensemble forecasts, developed at the US National Weather Service (NWS) and the Meteorological Service of Canada (MSC). When combined, these tools (a) provide weather forecast guidance for the 1-14 day period that is of higher quality than the currently available operational guidance based on either of the two sets of tools separately; and (b) make a set of forecasts that are seamless across the national boundaries over North America, between Mexico and the US, and between the US and Canada. As a first step in the development of the NAEFS system, the two ensemble generating centers, the National Centers for Environmental Prediction (NCEP) of NWS and the Canadian Meteorological Center (CMC) of MSC started exchanging their ensemble forecast data on the operational basis in September 2004. First NAEFS probabilistic products have been implemented at NCEP in February 2006. The enhanced weather forecast products are generated based on the joint ensemble which has been undergone a statistical post-processing to reduce their systematic errors.

5. Summary of 6th NAEFS workshop1-3 May, 2012 Monterey, CA 6th NAEFS workshop was held in Monterey, CA during 1-3 May 2012. There were about 50 scientists to attend this workshop whose are from Meteorological Service of Canada, Mexico Meteorological Service, UKMet, NAVY, AFWA and NOAA. • Following topics have been presented and discussed during workshop: • Review the current status of the contribution of each NWP center to NAEFS • For each NWP center, present plans for future model and product updates, for both the base models and ensemble system (including regional ensembles) • Decide on coordination of plans for the overall future NAEFS ensemble and products (added variables, data transfer for increased resolution grids, FNMOC ensemble added to NAEFS, especially for mesoscale ensemble-NAEFS-LAM) • Learn about current operational uses of ensemble forecast guidance, including military and civilian applications.

6. 7th NAEFS Workshop in Montreal, Canada • Time: 17-19 June 2014 • Locations: • 17-18 June – Biosphere, Montreal, Canada • 19 June – CMC, Dorval, Canada • Co-chairs: Andre Methot and Yuejian Zhu • Topics (or sessions) • Status and plan of Global ensemble forecast systems; • Operational data management and distribution; • Ensemble verification and validation metrics; • Reforecast, bias correction and post process; • Regional ensemble and data exchange; • Wave ensembles; • Integration of ensemble in forecasts: user feedback and recommendation; • Products – hazard weather, high impact weather and diagnostic variables; • Open discussion of the NAEFS research, development, implementation and operation plan

7. NAEFS Current Status Updated: November 18th 2014

8. Milestones • Implementations • First NAEFS implementation – bias correction Version 1.00 - May 30 2006 • NAEFS follow up implementation – CONUS downscaling Version 2.00 - December 4 2007 • Alaska implementation – Alaska downscaling Version 3.00 - December 7 2010 • Implementation for CONUS/Alaska expansion Version 4.00 - April 8 2014 • Implementation 2.5km/3km NDGD products for CONUS/Alaska Version 5.00 – August 2015 • Applications: • NCEP/GEFS and NAEFS – at NWS • CMC/GEFS and NAEFS – at MSC • FNMOC/GEFS – at NAVY • NCEP/SREF – at NWS • Publications (or references): • Cui, B., Z. Toth, Y. Zhu, and D. Hou, D. Unger, and S. Beauregard, 2004: “ The Trade-off in Bias Correction between Using the Latest Analysis/Modeling System with a Short, versus an Older System with a Long Archive” The First THORPEX International Science Symposium. December 6-10, 2004, Montréal, Canada, World Meteorological Organization, P281-284. • Zhu, Y., and B. Cui, 2006: “GFS bias correction” [Document is available online] • Zhu, Y., B. Cui, and Z. Toth, 2007: “December 2007 upgrade of the NCEP Global Ensemble Forecast System (NAEFS)” [Document is available online] • Cui, B., Z. Toth, Y. Zhu and D. Hou, 2012: "Bias Correction For Global Ensemble Forecast" Weather and Forecasting, Vol. 27 396-410 • Cui, B., Y. Zhu , Z. Toth and D. Hou, 2013: "Development of Statistical Post-processor for NAEFS"Weather and Forecasting (In process) • Zhu, Y., and B. Cui, 2007: “December 2007 upgrade of the NCEP Global Ensemble Forecast System (NAEFS)” [Document is available online] • Zhu, Y, and Y. Luo, 2014: “Precipitation Calibration Based on Frequency Matching Method (FMM)”. Weather and Forecasting (doi: http://dx.doi.org/10.1175/WAF-D-13-00049.1) • Glahn, B., 2013: “A Comparison of Two Methods of Bias Correcting MOS Temperature and Dewpoint Forecasts”MDL office note, 13-1 • Guan, H, B. Cui and Y. Zhu, 2014: “Improvement of Statistical Post-processing Using GEFS Reforecast Information”, Weather and Forecasting (in process)

9. NAEFS Statistical Post-Processing System • Bias correction: • Bias corrected NCEP/CMC GEFS and NCEP/GFS forecast (up to 180 hrs) • Combine bias corrected NCEP/GFS and NCEP/GEFS ensemble forecasts • Dual resolution ensemble approach for short lead time • NCEP/GFS has higher weights at short lead time • NAEFS products (global) and downstream applications • Combine NCEP/GEFS (20m) and CMC/GEFS (20m) • Produce Ensemble mean, spread, mode, 10% 50%(median) and 90% probability forecast at 1*1 degree resolution • Climate anomaly (percentile) forecasts • Wave ensemble forecast system • Hydrological ensemble forecast system • Statistical downscaling • Use RTMA as reference - NDGD resolution (5km/6km), CONUS and Alaska • Generate mean, mode, 10%, 50%(median) and 90% probability forecasts

10. NAEFS bias corrected variables Last upgrade: April 8th 2014 - (bias correction) 10

11. NAEFS downscaling parameters and products Last Upgrade: April 8 2014(NDGD resolution) 11 All downscaled products are generated from 1*1 degree bias corrected fcst. globally Products include ensemble mean, spread, 10%, 50%, 90% and mode

12. 2. NAEFS bias correction

13. NAEFS bias correction For any given forecast f, it could express as : Where a is truth (or replaced as best analysis), e is systematic error and ℇ is random error Therefore, systematic error (or accumulated bias) could be a time average difference of forecast and truth: In fact, an accumulated bias is changed from time to time which means e is not exactly systematic error, and ℇ is not exactly random error based on the time window for average.

14. NAEFS Bias Correction 1). Bias Estimation:We assume a bias (b) for each lead-time (t) (6-hour interval up to 384 hours), each grid point (i, j) is defined as the different of best analysis (a) and forecast (f) at the same valid time (t0) which is up on latest available analysis. 2). Decaying Average (or Kalman Filter method): Average bias will be updated by considering prior period bias and current bias by using decaying average (or Kalman Filter method ) with weight coefficient (w).

15. NAEFS Bias Correction 3). Decaying Weight: Through many experiments for different weights (w= 0.01, 0.02, 0.05, 0.1 and etc…), and different parameters, and different lead times, overall, w equals to 0.02 has been used for GEFS bias correction which is mainly using past 50-60 days information (see figure).

16. NAEFS Bias Correction 4). Bias corrected forecast: The new (or bias corrected) forecast (F) will be generated by applying decaying average bias (B) to current raw forecast (f) for each lead time, at each grid point, and each parameter. Simple Accumulated Bias Assumption: Forecast and analysis (or observation) is fully correlated

17. NAEFS Bias Correction 5). Performance: The performance is estimated by applying NAEFS bias correction method. The bias is calculated at each grid point for raw forecast (f) and bias corrected forecast (F), then using decaying average method (w=0.02) to get current average bias, taking absolute bias for each grid point, each lead time to generate domain average absolute error (bias) which smaller value is better (see figure: example for Northern Hemisphere 2 meter temperature, decaying average (w=0.2) about 2 months period ended by April 27, 2007).

18. 500hPa height: 120 hours forecast (ini: 2006043000) Shaded: left – raw bias right – bias after correction Forecast becomes worse after bias correction

19. 2 meter temperature: 120 hours forecast (ini: 2006043000) Shaded: left – raw bias right – bias after correction Positive bias Negative bias

20. Comparison of raw and bias corrected T2m for Summer and Fall 2011 RMS and Spread (NH – Fall) RMS and Spread (NH – Summer) Get worse from bias correction ME and ABSE (NH – Summer) ME and ABSE (NH – Fall) Raw forecast is bias free No carry on bias from summer

21. NAEFS Bias Correction • Several questions left behind • The correlation of prior joint sample • Assume samples are fully correlated • The weight to calculate decaying average • Optimum weights are functions of geographic and forecast lead times (should be) • Currently, w (weight) is fixed (w=0.02) • Systematic error for seasonal • Current method is lagged for seasonal information • 2nd moment adjustment for current method • Slightly adjust for CMC’s ensembles • N/A for single model ensemble

22. 3. Using ensemble reforecast

23. GEFS Reforecast Using Reforecast Data Configurations (Hamill et al, 2013) • Model version • GFS v9.01 – last implement – May 2011 • GEFS v9.0 – last implement – Feb. 2012 • Resolutions • Horizontal – T254 (0-192hrs– 55km) T190 (192-384hrs – 70km) • Vertical – L42 hybrid levels • Initial conditions • CFS reanalysis • ETR for initial perturbations • Memberships • 00UTC - 10 perturbations and 1 control • Output frequency and resolutions • Every 6-hrs, out to 16 days • Most variables with 1*1 degree • Data is available • 1985 - current • Bias over 24 years (24X1=24) • 25 years (25x1=25) • Bias over 25 years and 31day window • (25x31) • Bias over recent 2, 5, 10, and 25 years • within a window of 31day • (2x31, 5x31, 10x31, 25x31) • Bias over 25 years with a sample interval • of 7days within a window of 31days • and 61days (~25x4 and ~25x8) 1985 1986 . . . day-15 day day+15 2009 2010

24. Using 25-year reforecast bias (1985-2009) to calibrate 2010forecast Summer 2010 Solid line – RMSE Dash line - Spread Winter 2010 Spring 2010 Fall 2010

25. Using 24-year reforecast bias (1985-2008) to calibrate 2009forecast Spring 2009 Winter 2009 Decaying averages are not good except for day 1-2 2% decaying is best for all lead time Another year Summer 2009 Fall 2009 Decaying average is equal good as reforecast, except for week-2 forecast Decaying averages are not good except for day 1-3

26. Reforecast bias (2-m temperature) warm bias cold bias cold bias Perfect bias

27. T2m calibration for different reforecast sample sizes sensitivity on the number of training years (2, 5, 10, and 25 years) sensitivity on the interval of training sample (1 day and 7 days) Skill for 25y31d’s running mean is the best. 25y31d’s thinning mean (every 7 days) is very similar to 25y31d’s running mean. 25y31d’s thinning mean can be a candidate to reduce computational expense and keep a broader diversity of weather scenario!!! Long training period (10 or 25 years) is necessary to help avoid a large impact to bias correction from a extreme year case and keep a broader diversity of weather scenario!!

28. Using reforecast to improve current bias corrected product r could be estimated by linear regression from joint samples, the joint sample mean could be generated from decaying average (Kalman Filter average) for easy forward. Bias corrected forecast: The new (or bias corrected) forecast (F) will be generated by applying decaying average bias (B) and reforecast bias (b) to current raw forecast (f) for each lead time, at each grid point, and each parameter. Additional term (bias from reforecast) is added if r (correlation coefficient) is not equal one. This adjustment is expected to benefit for longer lead time forecast

29. Using reforecast to improve current bias corrected product (24-hr forecast, 2010 ) Spring Summer Winter Fall Perfect bias

30. Using reforecast to improve current bias corrected product (240-hr forecast, 2010 ) Spring Summer Winter Perfect bias Fall

31. Using 25-year reforecast bias (1985-2009) to calibrate 2010forecast 500hPa height Winter 2010 Very difficult to improve the skills Spring 2010

32. 4. 2nd moment adjustment

33. Make Up This Deficient ? Improving surface perturbations Or Post processing Bias correction does not change the spreads

34. 2nd moment adjustment 1st moment adjusted forecast 2nd moment adj. R=1 if E=0 Ensemble skill Estimated by decaying averaging Ensemble spread

35. After 2nd moment adjustment Winter Spring Almost Perfect!!! Summer Fall

36. Surface Temperature (T2m) for Winter (Dec. 2009 – Feb. 2010), 120-hr Forecast 25% spread increased SPREAD 7% RMSE reduced RMSE SPREAD/RMSE

37. Improving Under-dispersion for North America (Dec. 2009 – Nov. 2010) SPP2 RAW Winter Spring Summer Autumn RAW improving under-dispersion for all seasons with a maximum benefit in summer. Summer Winter Spring Autumn

38. 5. Multi-model ensemble

39. Bayesian Model Average Weights and standard deviations for each model (k - ensemble member) at step j Sum of (s,t) represents the numbers of obs. Finally, the BMA predictive variance is Between-forecast variance Within-forecast variance It is good for perfect bias corrected forecast, Or bias-free ensemble forecast, but we do not

40. Model 2 Model 1 Model 3 Courtesy of Dr. Veenhuis

41. Flow Chart of Recursive Bayesian Model Process (RBMP) Bias free Ens forecasts Observation or Best analysis Err and sprd Variance Weights Prior Err and sprd Prior weights Prior variance BMA, Decaying process, and adjustment New weights New variance New Err and sprd • (We thanks to Dr. Veenhuis for allowing us to adopt his BMA codes). 2nd moment adjustment Adjusted PDF

42. T2m for summer and Fall 2014

43. NUOPCIBC – simple combine three bias corrected ensembles DCBMA – decaying based Bayesian Model Average RBMP – Recursive Bayesian Model Process (built in decaying average and internal 2nd-moment adjustment) Solid line – RMS error Dash line - Spread Over-dispersion The results demonstrate: BMA could improve 3 ensemble’s mean, but spread could be over if original spread is larger RBMP could keep similar BMA average future, but 2nd moment will be adjusted internally All important time average quantities are decaying average (or recursive – save storage)

44. ROC are better for all leads Summer 2013

45. NUOPCIBC – simple combine three bias corrected ensembles DCBMA – decaying based Bayesian Model Average RBMP – Recursive Bayesian Model Process (built in decaying average and internal 2nd-moment adjustment) Solid line – RMS error Dash line - Spread Over-dispersion The results demonstrate: BMA could improve 3 ensemble’s mean, but spread could be over if original spread is larger RBMP could keep similar BMA average future, but 2nd moment will be adjusted internally All important time average quantities are decaying average (or recursive – save storage)

46. ROC are better for all leads Slightly degradation Fall 2013 Don’t understand this over-dispersion

47. U10m for summer and Fall 2014 Assume V10m has the similar statistics as U10m The difference of error/skills is very smaller, it is very hard to identify the improvement. Need to consider to have additional plots for difference of NUOPC and DCBMA and NUOPC and RBMA

48. RBMP is perfect for 1st moment and 2nd moment adjustment, especially for 2nd moment.

49. NH CRPS NH ROC Summer 2013 Tropical Southern hemisphere

50. The same conclusion as summer Perfect for 2nd moment adjustment