# Presentation: 416913

## Introduction of KMA statistic model and ensemble system

- Korea Meteorological Administration
- Numerical Weather Prediction Division
**Joo-Hyung Son**

**Notes:**
Hello, I’m Joo-Hyung Son from KMA
Today, I’d like to introduce the statistic model and ensemble model of KMA briefly. 

## Statistical models

**PPM (Perfect Prognostic Method)**
- Daily Max/Min and midnight temperature
- Probability of Precipitation
**MOS (Model Output Statistics)**
- Digital Forecast
**KF(Kalman Filtering)/DLM(Dynamic Linear Model)**
- Daily Max/Min Temperature
- 3 hourly temperature
- Daily Max/Min Temperature of 10 days

**Notes:**
At first, I’ll talk about the statistic model. Statistic models of KMA consist of 3 method, PPM,  MOS and KF, DLM. As you know, DLM is the just improved KF method. We use the PPM for daily Max/Min temperature and midnight temperature, and also generate probability of precipitation. 
Actually MOS is not developed in Numerical Weather prediction division which I’m working. This is developed for Digital Forecast in Digital Forecast division. They are generating wind, precipitation and other parameters with MOS as well as temperature. 
Kalman Filter is using just for temperature, daily max min temperature for short time, 3 hly temperature, daily max/min temperature for 10days.

## Statistical models

- KF

- PPM

- Max/Min Temp

- Max/Min Temp
- PoP

- GDAPS

- DLM

- DLM

- Max/Min Temp

- 3hr Temp

**Statistical models**

- PPM

- KF

- RDLM

- GDLM

**Notes:**
Each model is called just statistic method. So we call just PPM, KF, RDLM, and GDLM. RDLM mean the DLM method using regional model, so GDLM is global. 
PPM. KF, RDLM is the guidance for short range forecast. So the model for the implementation of statistic model is RDAPS, this is the regional model. Whereas the GDLM is for the week forecast. So it’s from GDAPS, global model. 

## PPM for Max/Min Temp

**Predictant**
- 00 UTC : +1(00UTC, Max/Min)
- 12UTC : +1(Max), +2(00UTC, Min)
**Forecast regions**
- 70 sites in Korea
**Model development**
- May 1, 1988 – Feb 28, 1992 (4 years)
- Regional reanalysis of JMA
- Climate data of 70 weather sites

**Notes:**
Let’s see the ppm.
We generate two times a day. 0 and 12 UTC for 70 weather stations in Korea with PPM. This model was developed 1998. The used data for generating was the regional reanalysis of Japan and climate data of each weather sites for 4 years from 1988 and 1992. When the model developed, we could not use MOS. Because we didn’t have enough model data.

## PPM Model structure for Max/Min Temp

**Temp (t) = A + B*obs(0) +   {C****i*****model predictor****i****(t)}**
**A, B, C****i**** (i=1,2,...,n): fixed coefficients**

**Forecast equation**

- 1000, 850, 700,500,400,300hPa
- Wind speed, direction, Temperature
- Dewpoint temp,
-  Height et al. from RDAPS
- Observation, climate

- Predictant

- Max/Min and 00LST temperature
- of 70 sights

- predictor

- Forecast eqs for
- each season, sights

- predictor

- predictant

**Notes:**
As you know, PPM is very simple method to use the multiple regression method. Forecast equation is like this. A, B and C are constant coefficients. 
Obs and model_I are the predictors.Forecast equations are separated according to the sights and season. Predictants are Max/Min and 0 LST temperature as I said. 

## PPM Predictors

-  select a group of predictors which explain predictant(temperature) well  from 44 predictors <method: forward-backward selection>
-  the number of the predictors of each seasonal and regional forecast equations are ranged from 5 to 10

| Main predictors for +12hr Max temp |  |  |  | Main predictors for +24hr Max temp |  |  |  |
| --- | --- | --- | --- | --- | --- | --- | --- |
| Spring | Summer | Fall | Winter | Spring | Summer | Fall | Winter |
| T850, CLMT,
OBS, VOR8
PCWT | OBS,
CTOP
VOR8
RH50,
PCWT | T85, CLMT,
OBS,
TTD8,
VORS | T85,
CLMT,
OBS,
TTDB,
VOR8 | T85,
CLMT,
PCWT,
VOR8,
OBS | OBS,
VORS,
T85,
CTOP,
RH70 | T85,
CLMT,
TTD8,
PCWT, VOR8 | T85,
CLMT,
VOR8,
TTD8,
70Q4 |
| Main predictors for +12hr Min temp |  |  |  | Main predictors for +24hr Min temp |  |  |  |
| Spring | Summer | Fall | Winter | Spring  | Summer | Fall | Winter |
| CLMT,
OBS,
T85,
KYID,
TTD8 | OBS,
T85,
CLMT,
VOR8,
TTD8 | CLMT,
OBS,
TTD8,
VOR8,
T85 | T85,
OBS,
TTD8,
TAD8,
VOR8 | CLMT,
T85,
KYID,
TTD8,
S70 | OBS,
T85,
CLMT,
VOR8,
PCWT | CLMT,
TTD8,
PCWT,
T85,
OBS | T85,
OBS,
TTD8,
VOR8,
TAD8 |

-  OBS: observation, CLMT: climate,  PCWT: virtual prediction, VOR: vorticity, TAD: temperature advection, KYID: KY index

**Notes:**
We selected predictors which explain predictant well from 44 predictor groups with forward-backward selection method. 
The number of the predictors of each seasonal and regional forecast equations are ranged from 5 to 10.
These are principle predictors for forecast
The most frequent factor is the 850hpa temperature and climate.
These are the principle predictors for each forecast. 

## PPM Model structure for PoP

- 1000, 850, 700,500,400,300hPa
- Wind speed, direction, Temperature
- Dewpoint temp,
-  Height et al. from RDAPS
- Observation, climate

- Predictant

- PoP of 18 regions

- predictor

- Forecast eqs for
- Each region according
- to warm and cold season

- predictor

- predictant

**Temp (t) = A + B*obs(0) +   {C****i*****model****i****(t)}**
**A, B, C****i**** (i=1,2,...,n): fixed coefficients**

**Forecast equation**

**Notes:**
The equation for PoP is same structure with temperature. But the equation was generated for 18 regions according to the warm and cold season. 

## PPM Predictors

-  PoP
-    the number of sites observed precipitation in the region
-               Total number of sites in the region
-  18 regions :
-        24 region by cluster analysis
-         (Moon(1990))
-        + forecast  experiment

- 18 regions for forecast of PoP

-  the forecast equations are developed according to the warm(April-September) and cold(October-March) season and each regions.

| principle predictors for PoP |  |
| --- | --- |
| Warm season | DWL, VR850, QA700, VV700, S850, RH500 |
| Cold season | DWL, VV850, RH850, VR850, S850, 7Q4 |

**Notes:**
We use the definition of PoP, the ratio of the number of sites observed precipitation in the region and the total number of sites in the same region.
The forecast region was divided by cluster analysis and forecast experiment. 
Forecast equations are developed according to the arm and cold season and each regions. 
Next is the principle predictors for PoP
(DWL: depth of wet layer)

## Slide 9

**Notes:**
These are the graphic we are displaying in web-site.

## KF for Max/Min Temp

**Predictant**
- 00 UTC : +1(Min/Max), +2(Min)
- 12 UTC : +1(Max), +2(Min/Max)
**Forecast regions**
- 40 in Korea,
- 32 in North Korea, China, Japan

**Notes:**
Klamn Filter is also generated two times a day. This is forecast 40 Korea weather sites and 32 principle cities of North Korea, China, and Japan.

## Kalman Filter algorithm

- KF for Max/Min Temp

- vt~N(0,Vt): observation noise

- wt~N(0,Wt) : process noise

-                     1
- Ft =           RDAPS
-          Latest Obs temp

- Gt = 1

- V0 = 2

-             4/365     0       0
-   W0=      0    1/365    0
-               0        0    1/365

**Notes:**
As you know, Kalman filter is the iteration between observation and system equation. 
In our system, we are fixing the initial observation noise and process noise covariance and Gt. Ft consists of 1, RDAPS model forecast and latest observation temperature. 

Wt : process noise coveriance
Vt: observation noise covariance

## Slide 12

**Notes:**
This is the KF graphic.

## DLM(Dynamic Linear Model)

- DLM
- Improved Kalman Filter algorithm
- Weights(regression coefficient) are modified according to the prior condition with time.

**Notes:**
DLM is the improved kalman filter algorithm. Liner regression has the fixed coefficients, whereas DLM has dynamic coefficients according to the prior condition with time.  

## DLM

- vt~N(0,Vt)

- wt~N(0,Wt)

- Use the updating algorithm to estimate Wt with time
- Find appropriate Wt increasing discount factor(0<delta<1) from 0.01 to 1 with interval 0.01
- the discount factor is selected when RMSE between observation and forecast is the lowest

- DLM(Dynamic Linear Model)

**Notes:**
In the DLM, Wt is changed. We use the updating algorithm to estimate Wt with time.
To find the appropriate Wt increase discount factor from 0.01 to 1 with interval 0.01
Discount factor is selected when the REMS between observation and forecast is the lowest.
Wt : process noise coveriance
Vt: observation noise covariance

## RDLM(Regional DLM)
3hourly forecast up to 48hr
RDAPS 
38 sites

GDLM(Global DLM)
Max/Min temp for 10 days
GDAPS
38 sites

- 3hourly forecast up to 48hr
- RDAPS
- 38 sites
**GDLM(Global DLM)**
- Max/Min temp for 10 days
- GDAPS
- 38 sites

- DLM(Dynamic Linear Model)

**Notes:**
Two kind of forecasts are generated with DLM method. RDLM is 3houly forecast up to 48hr for 38 sites.
GDLM is the Max/Min temperature for 10days at same sites. 

## Slide 16

**Notes:**
This is the RDLM forecast. It just shows the time series. 

## 
Ensemble Prediction System

- Ensemble Prediction System

**Notes:**
Next, I’ll talk about ensemble forecast system briefly. Actually really briefly.
I was not responsible for the ensemble forecast system while I was in Korea. So I don’t know well about this.

##    KMA Ensemble Prediction System 

|  | GBEPS 1.1.1
~ GBEPS 1.2.1 | GBEPS 2.1.1
~ GBEPS 2.3.1 | GBEPS 2.3.1.1
~ GBEPS 2.3.1.2  |
| --- | --- | --- | --- |
| Operation period | 2001.3.1∼2003.10.31 | From 2003.11.1 | From 2005.2. |
| Data assimilation | 2dOI → 3dOI | 3dOI → 3dVar | 3dVar |
| Model | GDAPS T106L21 | GDAPS T106L30 | GDAPS T106L30 |
| Vertical resolution | 21 levels | 30 levels | 30 levels |
| Perturbation method | Breeding | Breeding → Breeding + Factor Rotation | Breeding → Breeding + Factor Rotation |
| Target area (BV) | Global | Northern Hemisphere | Northern Hemisphere |
| Lead time | 10 days | 8 days | 8 days |
| Ensemble members | 17 (16 members + 1 control) |  | 17 members |

**Notes:**
Our ensemble prediction system based on the breeding method with the global spectral model T106L30, has been on operation since March 2001.

After porting 2nd super computer, we run two times a day since the computer resource is improved. So we can generate 34 members to make several ensemble products.

## Schematic diagram

**Schematic diagram**

-   The global spectral model T106L30 with the slightly different initial conditions run 17 times.
-   Both perturbed analysis and control analysis are projected to 24hours with the model, and departures from the control analysis at +24hours are scaled down to the norm of initial perturbations

## Schematic diagram

- D day + 12hr

- D+1 day

- D+1 day +12hr

**Schematic diagram**

-   17members could be similar each other because they are generated from the identical model, so this is to make different perturbation among the members manually.
-   In the new system, the factor rotation was added every alternative step.

**Notes:**
I don’t know well the perturbation method. What I can explain is that we added the factor rotation in the middle of breeding method. 
The purpose of this is to generate perturbation which can distinguish each member more apparently

멤버를 17개 만들어도 같은 모델에서 성장시키기 때문에 서로 유사한 패턴을 가지기 쉽기 때문에 12시간   브리딩 한 멤버들을
orthgonal 하게 (서로 직교성을 가지게)  회전시켜 인위적으로 서로 다른 섭동장을 만들어서  다시 12시간 브리딩 시킨답니다. 

## 2005. 6. 11

**2005. 6. 11**

**old(cray-before)**

**NEW (cray-frot)**

**Notes:**
There are some results to show the improvement of the ensemble. 
The spread of short lead time is extended compare with old system.

## 2005. 6. 11

**2005. 6. 11**

**old(cray-before)**

**NEW (cray-frot)**

**Notes:**
And according this figure, it seems that the RMSE is also improved. Compare the 5 days RMSE 

##   EPS products (http://190.1.20.56)

**Notes:**
We are offering some products in web-sites. ~

**Mean and Spread**

**Notes:**
This is the ensemble mean and spread 500 hPa hight, MSLP, 850hPa temperature, respectively.

##  Spaghetti

** ****Spaghetti**

**Notes:**
Spaghetti diagram can plot just one contour for given 500hPa height from each member on the same map. 
If we loop these plot by forecast hour, we can discriminate the uncertainty. We can immediately discern where the members tend to lie close together is more predictability. 

 all the members lying orderly in the early  forecast hour
 and then  slowly diverging from each other as lead time increases.

##  Spaghetti ( with ensemble spread)

- 5640m

- 5520m

##  Stamp map  

- • display the global model, mean and standard deviation and spaghetti as well as each member.

**Notes:**
앙상블 예보는 평균장과 편차장을 제공하는 것도 중요하지만 모든 멤버의 시나리오를 함께 표출하여 예보에 활용하는 것도 중요하다. 이에 2004년부터 500hPa 고도에 대하여 17개 앙상블 예보장을 우표처럼 작게 붙여 한 지면에 표출하는 Stamp map을 수치예보과 인트라넷에 제공하여왔다.   앙상블 모든 멤버들에 대한 시나리오 이외에 고해상도 전지구예보(GDAPS  T213)를 함께 표출함으로써 extreme 시나리오의 예측 가능성을 검토 가능하게 하였고, 특히 stamp map의 효과적인 활용을 위해 앙상블 평균장, 편차장 및 스파게티도 참고 자료로 표출하도록 개선되었다. 또한 5520m, 5640m, 5760m의 세 고도를 선택하여 색깔을 달리 표현함으로써 시나리오간의 각기 다른 기압 패턴을 파악하기에 용이하게 하였다. 
한편 500hPa 고도와 함께 해면기압에 대한 stamp map도 표출하고 일기도의 영역과 크기를 개선할 예정이다. 

##  12-hour precipitation  >  given  thresholds 
: 1, 5, 10mm for winter season
: 1, 10, 50mm for other seasons 
The probability 

These probability maps are used for the early warning guidance of severe weather.  

- : 1, 5, 10mm for winter season
- : 1, 10, 50mm for other seasons
- The probability
- These probability maps are used for the early warning guidance of severe weather.

-   Categorical PoP

**Notes:**
Categorical PoP is calculated from 12hour precipitation given thresholds. Thresholds are 1, 5, and 10mm for winter, otherwise 1, 10, 50mm
The formula is very simple. We define the probability is the ratio the number of member above threshold and the total number of ensemble member. 

##   Categorical PoP

-   Categorical PoP

**Notes:**
These are showing the categorical PoP. 

##  Surface maximum wind  >  10m/s, 14m/s
 The probability 

These probability maps are used for the early warning guidance of severe weather. 

-  The probability
- These probability maps are used for the early warning guidance of severe weather.

-   Probability of Surface Max Wind

**Notes:**
The probability of surface max wind is the same with categorical PoP. The thresholds are 10 and 14m/s.

##   Probability of Surface Max Wind

-   Probability of Surface Max Wind

**Notes:**

This is the figure of probability of Max Wind

##     Time series of Probability 

-  **Precipitation**
-    12hr accumul >= 1mm
-    12hr accumul >= 10mm
-    12hr accumul >= 50mm
-   **Surface Max Wind**
-     sfc wind  >=  10m/s
-     sfc wind  >=  14m/s
-  **Principle cities**
-     Seoul, Daegu, Daejeon
-     Busan et al.

**Notes:**
We are also plot the probability with time series for principle cities.

##    EPSgram

**  ****Time series of primary cities**

**Notes:**
Precipitation, 10m wind, daily max wind

## Ensemble Plumes

**Ensemble Plumes**

**Time series of 8-day forecast at cities**
**The dispersion of members with forecast evolution**** **
**Variable : Pmsl, 500H, 850 T**

## Hwangsa (yellow sand) trajectory 

## Typhoon Strike Probability Map by EPS

- Thank you

##    Factor analysis 

- Factor analysis
- Factor analysis is a statistical technique to explain the most of the variability among a number of observable random variables in terms of a smaller number of unobservable random variables called factors
- Factor rotation
- Factor rotation is to find a parameterization in which each variable has only a small number of large loadings. That is, each variable is affected by a small number of factors, preferably only one. This can often make it easier to interpret what the factors represent.