Thursday, May 19, 2016

Application of the Markov Process Rainfall Model

Posted by Jason Lillywhite

If you have visited our Model Library lately, you might have noticed that we have a nice little example model that demonstrates the use of a Markov process to simulate daily rainfall. You need to specify some key statistical inputs that have some basis on historic data. How do you develop these inputs? How do you know if the Markov model is realistic? I thought it would be helpful to show how this simple example might be applied in a real-world project and try to answer those questions.

*Note that I made changes to the results on 5/20/2016 after I used GoldSim's optimization function to better calibrate the rate variability.

Review the Historic Record

The first step of the process is to obtain daily precipitation records for your project site location. Make sure you have at least 10 years of daily data. Much less than this will make it quite difficult to validate and calibrate your model. For my example, I obtained daily rainfall data from Mt. Hotham, Victoria, Australia for the period of Jan 2002 thru Jan 2014. Using nothing more than a single time series element in a GoldSim model, I'm able to use reporting periods and time series shifting to quickly generate a summary of the historic data that provides a good overview of the data. Below is a screen capture showing how I've set up this simple model:

Below are 3 different views of the precipitation time history after running the model for the simulation duration.
View of rainfall time history

View of rainfall time history using monthly reporting periods

View of rainfall time history using annual reporting periods

By changing the simulation settings to run for a single year, Monte Carlo for 12 realizations where each realization represents each year from 2002-2014, we can enable the time shifting option of the time series and then quickly calculate the probability distribution of annual cumulative rainfall.
Probability Distribution Plot of Annual Rainfall Totals

The Markov Process GoldSim Model

Now that understand the basic characteristics of our rainfall record, let's take a look at the input requirements of the Markov process example model. If you are following along, please go to our Model Library and download the model from here. A screen capture of the model is shown below. Basically, this model works by switching between dry and wet states during the course of the simulation. It cannot be in wet and dry at the same time and while in one state, the other is inactive. As soon as one becomes active, the other is deactivated. The triggering of wet/dry states is driven by the probability of it being wet at any given time and the average duration of the "wet" state or period in which it is raining. The states are controlled in GoldSim using conditional containers, which allow you to control their activity dynamically using time-based triggers.

The 3 elements on the left are primarily what are used to drive this model. The "WetLength" is a scalar constant value that specifies the mean duration of wet periods. The next input is a table of monthly average rain depths and the final table is monthly average probability of a wet state occurring. This information can be obtained from an analysis of the historic record. 

Analyze the Historic Record

I used another GoldSim model to analyze the historic record but it is not the only way to do it. I could have used Excel, Matlab, or any other time series analysis tool to do this. If you are interested in obtaining a copy of the model I used to do the analysis, please contact us. To calculate the average wet duration, I developed triggering logic that uses criteria to determine the start and end time of each rainy period then I capture those duration values in an array, which I then use to calculate the overall mean duration of all the storms that occurred in the historic record. The trigger of a beginning of a storm requires that it was dry and that the amount of cumulative rain since rain began exceeds a threshold value. The trigger to mark the end of a storm requires that the rain has dropped below a threshold and held that value for at least 1 day. I calculated the average storm duration for each month. Below is the table showing these values.

I used these values for the "WetLength" input data element. To capture the changes for each month, I refer to the data using the vector index lookup: "WetLength[month]".

The "PrecipHist" table is the monthly average rain depth from the historic record. I obtained this from GoldSim using an Integrator with the "type" set to a vector of Months and the rate of change input set to: vector(if(row = Month, Rainfall_Data, 0 mm/hr)). This will sum the rain for each month as the simulation walks through all the years on record. Then I just averaged those values over the simulation period to obtain the numbers for the Markov model.

Finally, I used a Selector element to determine whether it was raining or not and then calculated the mean probability of a wet day for each month.

I also used a submodel in GoldSim to determine the day to day variability of the precipitation rate by plotting a probability distribution of each value at all the time points. I found that this distribution fits fairly well with a Log Pearson III distribution with location, scale, and shape factors of 0.03, 1.4, and 0.66 respectively. I determined these values using GoldSim's built-in optimization routine. The objective function was to minimize the difference between the historic and simulated rainfall depth and peak rate.

Results Validation

Using this information, I configured the Markov process model and also added to it the historic data for validation and comparison. I'm running the Markov model for 1 calendar year, 100 realizations. After inspecting a few realizations of daily precipitation time series, I found the model to be producing fairly reasonable results as can be seen in the plots below, which show a single realization and the mean of all realizations for daily values.

It should be noted that for the mean time history, you would expect the historic values to look more volatile because it is drawing from far fewer possibilities (only 12 for years 2002-2014) while the Markov model is running with 100 possibilities. This will cause the mean time history to show less change from day to day as more possibilities are taken into account. So the Markov model appears to produce sensible results based on the daily time histories but what about the peak, cumulative annual, and monthly statistics? Below is the plot showing the probability distribution of annual peak rainfall rates. The range appears to be inline but the Markov model is showing much higher peak rain rates. Since the lower values appear to match well, .
I truncated the Log Pearson III distribution at a peak scale factor of 10, which appears to work well and still allows for a good fit on the total annual rainfall. See the plot below, which shows the peak rain rates after I made this adjustment.

The match between the simulated and historic values for total rain appears to be quite good. 

The monthly mean, 10%, and 90% statistics do appear to match quite nicely. This is to be expected because the model is being driven by the monthly mean precipitation values for each month. This also means the model is limited by the few years of record we are given to do this analysis.


The really nice thing about this model is that we can pretty quickly generate stochastic rainfall inputs to our water management models based on historic data without too much effort. By going through the type of analysis shown above, you can provide a sense of confidence to your audience and the end-users of your models while also adding robust capability to your model. Not only will you incorporate the historic possibilities but also those not yet recorded and therefore provide a way to visualize the otherwise unexpected. That is, after all, what GoldSim is for - to help you make decisions in an uncertain world. 

I hope you find this article helpful. Please let us know if you have any questions by submitting a request in our Help Center or writing in our Forum.


No comments:

Post a Comment