From: Kirk on
I would like to double check and that I am using the 'mvnrnd' command properly, but also I am trying to better understand the details of the 'mvnrnd' command. I will try and be as concise as I can while still providing enough detail to show what I am doing.

I am using 'mvnrnd' to generate generate random climate arrays of 4 monthly climate variables (max temp, min temp, solar radiation, and precipitation). These arrays in turn, drive a ecological model. Obviously these are highly correlated variables.

First I load a 12 month x 4 variable marix of means (µ), and 4 variable x 4 variable x 12 month matrix of covariances (sigma). These are derived from 100 years of measured monthly climate station data. I use 'repmat' to "grow" the µ and sigma matrices to the desired length (let's say 1000 years). I also using a simple ramping function on the 2 temperature means to simulate warming on the final 200 years of the array of repeated means. I am simulating increased variability of the climate by using a multiplier to ramp the final 200 years of the covariance matrix. This multiplier acts on all elements of the covariance matrix thereby maintaining the covariance relationships among tmax, tmin, rad, and precip. In the end of have a array of µ and an array of sigmas with some combination of climate and increased variability on the final 200 years of the array. I then run 'mvnrnd' on the full
1000 years of repeated µ and sigma to yield a 1000 year array of random climate variables.

My questions and concerns all circle around the idea of how 'mvnrnd' maintains the covariance relationships its input parameters.

1. I understand that 'mvnrnd' maintains the covariances across variables within a given month, but, does 'mvnrnd' also maintain covariances from month to month along the same variable? To say that another way... it is obvious that the variables tmax, tmin, radiation, and precipitation will covary at any time-step (which is to say, a dry July is likely to be a sunny July), However, it is also likely that a hotter than average July is more likely to be followed by a hotter than average August. Hence my question, Does mvnrnd take this temporal covariance into account?

2. Should I be applying the the variance multiplier to all elements of the covariance matrix, or only to the diagonal (variances of the covariance matirx)? or do I missunderstand the nature of the elements of a covariance marix?

Any thoughts or ideas that addresses these concerns are much appreciated.
From: Peter Perkins on
Kirk wrote:

> 1. I understand that 'mvnrnd' maintains the covariances across variables
> within a given month, but, does 'mvnrnd' also maintain covariances from
> month to month along the same variable?

If I understand your description correctly, each month's value will be a 1x4 vector. As you say, the covariance among those four variables will be the 4x4 cov matrix you specify (and the sample cov will approximate that), but the samples are INDEPENDENT from month to month. So, as you appear to be using MVNRND, the answer is no.

> To say that another way... it is
> obvious that the variables tmax, tmin, radiation, and precipitation will
> covary at any time-step (which is to say, a dry July is likely to be a
> sunny July), However, it is also likely that a hotter than average July
> is more likely to be followed by a hotter than average August.

In the absence of a trend in the mean, you'll need a more complicated cov structure to model that kind of thing. A moving average, perhaps.


> 2. Should I be applying the the variance multiplier to all elements of
> the covariance matrix, or only to the diagonal (variances of the
> covariance matirx)?

It sounds like what you want is to inflate the variances, while preserving the correlation. That implies that you need to scale the off-diag terms in the cov matrix. I believe what you should do is multiply the original cov matrix elementwise by sqrt(varscale)*sqrt(varscale'), assuming varscale is a col vector. If you have a scalar scale factor, then it's much simpler -- all you need do is scale the entire cov matrix.

Hope this helps.
From: Kirk on
Thank you for the reply Peter...

Peter Perkins <Peter.Perkins(a)MathRemoveThisWorks.com> wrote in message <hjlont$12j$1(a)fred.mathworks.com>...
> Kirk wrote:
>
> > 1. I understand that 'mvnrnd' maintains the covariances across variables
> > within a given month, but, does 'mvnrnd' also maintain covariances from
> > month to month along the same variable?
>
> If I understand your description correctly, each month's value will be a 1x4 vector. As you say, the covariance among those four variables will be the 4x4 cov matrix you specify (and the sample cov will approximate that), but the samples are INDEPENDENT from month to month. So, as you appear to be using MVNRND, the answer is no.
>

You understand correctly. And I suspected as much...

> > To say that another way... it is
> > obvious that the variables tmax, tmin, radiation, and precipitation will
> > covary at any time-step (which is to say, a dry July is likely to be a
> > sunny July), However, it is also likely that a hotter than average July
> > is more likely to be followed by a hotter than average August.
>
> In the absence of a trend in the mean, you'll need a more complicated cov structure to model that kind of thing. A moving average, perhaps.

I'm not sure how to handle this... As you imply, things get complicated. Perhaps my first approach should be to try and come up with a test for temporal autocorrelation. Although even this gets tricky, because I am sure that the periodicity of monthly temperatures in the measured climate station data will show a high degree of autocorrelation over the year. Perhaps being able to show a similar degree of autocorrelation in the mvnrnd modeled data would reduce some of my concerns?

>
>
> > 2. Should I be applying the the variance multiplier to all elements of
> > the covariance matrix, or only to the diagonal (variances of the
> > covariance matirx)?
>
> It sounds like what you want is to inflate the variances, while preserving the >correlation.

That is correct. Although to be more specific, I am gradually inflating the variances with a 200 year (2400 month) "ramp".

>That implies that you need to scale the off-diag terms in the cov matrix. I
>believe what you should do is multiply the original cov matrix elementwise by >sqrt(varscale)*sqrt(varscale'), assuming varscale is a col vector. If you have a scalar >scale factor, then it's much simpler -- all you need do is scale the entire cov matrix.

I think I have a scalar-scalar factor, but I'm not entirely sure what you mean here. Below is an example of my covariance matrix and an example multiplier matrix:

The covariance matrix:

climate_sigma <4x4x12 double>

val(:,:,1) =

1.0000 0.9362 0.9162 -0.1859
0.9362 1.0000 0.7189 -0.1899
0.9162 0.7189 1.0000 -0.1995
-0.1859 -0.1899 -0.1995 1.0000

val(:,:,2) =

1.0000 0.9391 0.8761 -0.2509
0.9391 1.0000 0.6589 -0.1961
0.8761 0.6589 1.0000 -0.3211
-0.2509 -0.1961 -0.3211 1.0000

..
..
..

val(:,:,12) =

1.0000 0.9061 0.6831 -0.2267
0.9061 1.0000 0.3790 -0.1647
0.6831 0.3790 1.0000 -0.2576
-0.2267 -0.1647 -0.2576 1.0000

(climate_sigma gets duplicated with "repmat" to become a <4x4x2400 double> named sigmaPost)

Here is the covariance multiplier:

covarInc7 <4x4x2400 double>

val(:,:,1) =

1.0025 1.0025 1.0025 1.0025
1.0025 1.0025 1.0025 1.0025
1.0025 1.0025 1.0025 1.0025
1.0025 1.0025 1.0025 1.0025


val(:,:,2) =

1.0050 1.0050 1.0050 1.0050
1.0050 1.0050 1.0050 1.0050
1.0050 1.0050 1.0050 1.0050
1.0050 1.0050 1.0050 1.0050

..
..
..

val(:,:,2400) =

7 7 7 7
7 7 7 7
7 7 7 7
7 7 7 7

Then I multiply to two matrices with the .* operator

sigmaPost.*covarInc7; % apply covriance increase here

I think this is what you are referring to when you use the phrase scalar-scalar multiplication?

>
> Hope this helps.
From: Peter Perkins on
On 1/26/2010 3:14 PM, Kirk wrote:

>> That implies that you need to scale the off-diag terms in the cov
>> matrix. I believe what you should do is multiply the original cov
>> matrix elementwise by >sqrt(varscale)*sqrt(varscale'), assuming
>> varscale is a col vector. If you have a scalar >scale factor, then
>> it's much simpler -- all you need do is scale the entire cov matrix.
>
> I think I have a scalar-scalar factor, but I'm not entirely sure what
> you mean here.

All I meant was that if you have a cov matrix S for X1, X2, and X3, then a cov matrix for a1*X1, a2*X2, and a3*X3 is ([a1;a2;a3]*[a1,a2,a3]).*S, and the variances along the diagonal of that are inflated by [a1^2,a2^2,a3^3]. The correlation matrix among a1*X1, a2*X2, and a3*X3 is the same as that for X1, X2, and X3.

And if a1==a2==a3, i.e., you have a scalar scale factor a, then all you need is a^2*S.

It's up to you what you mean by "variance multiplier", but I would think you would want to do something like the above.