Steve Markstrom: Hello!
This is Steve Markstrom as part of the PRMS training series, I am going to talk a little
bit about the parameter sensitivity program FAST and how we apply it to PRMS.
To start out with, FAST stands for, Fourier Amplitude Sensitivity Test.
It's a variance based global sensitivity algorithm for estimating first order partial variance
explained by each parameter and we're using the R software package.
The idea is that if amplitude of this sine curve represents range of the parameter values
that we want to test, we can vary the values of different parameters at different frequencies
over the trial space.
Along the X axis is some numbers of trials that we have to test.
A point along the X-axis represents a trial, for instance this is the value of P3 to use,
this is the value of P2 to use and this is the value of P1 to use for that particular
trial.
FAST is sampling the parameter space like that.
You run all those different trials and you compute an objective function for each trial
and you get some values that can be plotted like this and then you run a Fourier analysis
on it.
You can see the frequencies associated with those different parameter sets in the output.
This would be the power spectrum associated with the parameter P1, power spectrum with
parameter 2 and so on.
By looking at these different power spectrums, you can assign sensitivities to the parameters.
I have put together a little package called PRMS FAST.
This is not something I distribute, not something I have documented, this is just software I
use and if at the end of this presentation you feel like you are interested and capable
and willing to take on something like this on let me know.
Send us a request through our help system and perhaps I will share this with you.
It sits on top of a PRMS workspace.
I'll show you how that works.
Now let's run a little demo.
This is the standard ACF PRMS workspace but it also has the PRMS FAST stuff in here.
If we look at this ‘README_prmsFast’, bring that up, here are the steps.
First thing you have to do is install R and install the FAST package.
Hopefully you paid attention to Lauren’s LUCA presentation where she told you to make
a table like this, these are going to be some of the parameters that we use to test the
sensitivity on.
Then if you look at your workspace here, you can find this CSV file that you can open in
excel, that’s pretty easy.
What you do is in this first column you put the name of the parameter.
That’s the PRMS parameter.
You can look them up on user's manual and they're in there and then what you do is set
the range.
The minimum and the maximum and this is setting the amplitude of how FAST is going to vary
these parameters, so I could add more parameters or delete parameters.
But this is a fairly small set and it's nice for a demo.
That's the input.
The first program you're going to run is called ‘fast.R’.
Here it is and now I run it.
It’s reading that ‘sens_param.csv’ that we just looked at in excel.
Now what am going to do is come over here and here's ‘trials.csv’, so I just ran
that R program and it wrote this file right here.
The first column is the trial number and these are the values of the parameters that it wants
me to run.
FAST want me to make all these runs, it wants me to make 91 runs and set the parameters
to these values.
That's what I am going to do.
If I look at the README, it says run Java program ‘prmsFast’.
What ‘prmsFast’ does is it looks at this ‘trials.csv’ file and it writes a parameter
file and a control file for every one of those trials, it creates the trials.
I am just going to run that.
You can see that's what it's doing.
It wrote the control files and it's done and now.
If I go into the control directory, there is a directory ‘sens’, and you can see
I have 91 control files and if I go to the ‘input’ directory, you can see the trial
parameter files.
That Java program I just ran wrote all these files.
It also wrote this big batch file that makes all the runs.
If I fire that off it's going to start running all those trials.
There is also a Shell script for Linux and a SLURM script for our parallel cluster.
You can see it's starting to run and it's writing in the output directory.
We'll take a break here and I’ll wait for this to finish.
OK, the runs are done now and if we go and look at our README file, you can see that
the next step, the next step, is to make sure that the runs actually finished running.
And that’s what I’m going to do next.
If we go to our output directory, there you can see the runs that were just made by the
batch file that we ran.
Then if we look at the README again, the next step is run this R script ‘prmsSens.R’,
so I’ll start that now.
All right, so here’s the R file, ‘prmsSens’.
If we look down inside of this file, here you can see where we’re looping over each
of the trials reading the CSV file in the output directory.
And right here, this line right here, is where the objective function is computed.
This is a root mean squared error function.
If you wanted a different objective function, you could write your own function right here,
but this is what I am going to run.
So I am going to highlight this entire file and I’m going to run it.
It's reading the output files from the 91 trials.
Its computing the objective functions and then its running a Fourier frequency analysis
on that trial space.
When it's done we will look at the output.
OK, there it's done.
If we look in the main ACF directory, there are two csv files.
The objective function root means square R values, that’s this one here, but this is
what we really want to look at is the second one, ‘sens.csv’.
The first column lists the different output variables and the rows across the top lists
the sensitivity parameters.
It doesn't necessarily compute a sensitivity value for all of them.
If there is no delta then you're not going to get a sensitivity value.
That’s what these ‘NAs’ are.
This is the way to look at it.
Choose a variable, ‘basin_cfs’ for instance.
Then look across the corresponding row for the largest number, then go up that column
and see the parameter that influences this variable the most.
In this case, it’s ‘soilmoist_max’.
That was a little demonstration of how to run it for yourself.
Now, I’m going to close all of this stuff up.
I happen to run this on 35 parameters there they are for eight different output variables
there they are.
For seven different statistics, here are those.
For each of the hydrologic response unit of the national model, and there are those.
For ten years.
Here’s some output of that study.
You can see different objective functions across the top, autoregressive lag-1, coefficient
of variation and mean.
There are different hydrologic processes listed down the side, baseflow, ET, runoff, infiltration,
snowmelt, soil moisture, surface runoff, and interflow.
One thing that kind of struck me was how different they all seem to look.
Depending on what you're actually interested in and what's your objective function is you're
going to get some very different answers.
I don't have the scale on here but it goes from white to pink to green to dark blue.
White is a zero value and dark blue is a value of one because sensitivity runs from zero
to one.
I think Lauren showed this plot and this is kind of interesting.
For every HRU if you rank the processes by the most dominant to the least dominant this
plot on the top is showing the most dominant process, which for most of the countries seems
to be ET.
All of these processes are dominant somewhere which is a very interesting outcome and not
necessarily the way it had to be.
Then we did the same thing down here with the least dominant processes.
The idea would be is if you are going to calibrate in these areas you will want to think about
these processes and you could simplify your problem by removing these processes from your
either calibration or evaluation whatever you are doing.
That concludes this presentation.