从Matlab中相关的边际统一中提取:copularnd的替代方案? [英] Drawing from correlated marginal uniform in Matlab: Alternatives to copularnd?

查看:223
本文介绍了从Matlab中相关的边际统一中提取:copularnd的替代方案?的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

考虑两个随机变量XY均均匀地分布在[0,1]中并与相关性rho相关.

我想从Matlab中的联合分布中提取P实现.

A = copularnd('gaussian',rho,P);

唯一的方法吗?

解决方案

科普拉斯提供一种非常方便的方式对关节分布进行建模.仅说要X~U[0,1]Y~U[0,1]corr(X,Y)=rho不足以定义这两个随机变量之间的关系.只需通过替换不同的copula,您就可以构建都满足此条件的不同模型(最终选择最适合您的用例的模型).

除了了解系动词的基本知识外,您还需要了解各种相关类型,例如 linear ( Spearman / Kendall ),以及他们如何彼此相关.特别是,与线性相关不同,秩相关在应用单调变换时会保留.这是关键属性,它使我们能够轻松地将所需的均匀边际分布线性相关性转换为双变量正态线性相关性(或在copula中使用的其他类型的分布),这将是与copularnd的输入相关性. /p>

关于交叉验证,有一个很好的答案,准确描述了如何在您的情况下转换相关性.您还应该在使用等级相关系数上阅读此Matlab指南. 和copulas.

简而言之,要对所需的边际线性相关性进行建模,您需要将其转换为某种等级相关性(使用Spearman很方便,因为对于均匀分布,它等于Pearason相关性),这与您的法线相同(因为这是一个单调的转换).您所需要做的就是将法线的Spearman相关性转换为线性Person相关性,这将成为copularnd的输入.

希望该代码可以使其易于理解.显然,您不需要这些临时变量,但是应该使逻辑清晰:

rho_uni_Pearson = 0.7;
rho_uni_Spearman = rho_uni_Pearson;
rho_normal_Spearman = rho_uni_Spearman;
rho_normal_Pearson = 2*sin(rho_normal_Spearman*pi/6);

X = copularnd('gaussian', rho_normal_Pearson, 1e7);

得到的线性相关正是我们想要的(因为我们已经生成了一个非常大的样本):

corr(X(:,1), X(:,2));
ans =
    0.7000

请注意,对于二元正态copula,线性和秩相关之间的关系很容易,但是如果您使用其他copula,则关系可能会更加复杂. Matlab具有功能copulaparam,使您可以将秩相关转换为线性相关.因此,我们不用写上面的显式公式,而可以使用:

rho_normal_Pearson = copulaparam('Gaussian', rho_normal_Spearman, 'type', 'spearman')

现在我们已经学习了基础知识,让我们继续使用具有5自由度的t copula代替高斯copula:

nu = 5; % number of degrees of freedom of t distribution
rho_t_Pearson = copulaparam('t', 0.7, nu, 'type', 'spearman');

X = copularnd('t', rho_t_Pearson, nu, 1e7);

确保结果线性相关是我们想要的:

corr(X(:,1), X(:,2));
ans =
    0.6996

很容易看到,根据您选择的系动词,即使它们都为您提供了相同的线性相关性,所得的分布也可能会显着不同.研究人员将由您决定哪种模型最适合您的特定数据和问题.例如:

N = 500;
rho_Pearson = copulaparam('gaussian', 0.1, 'type', 'spearman');
X1 = copularnd('gaussian', rho_Pearson, N);
figure(); scatterhist(X1(:,1),X1(:,2)); title('gaussian');

nu = 1; % number of degrees of freedom of t distribution
rho_Pearson = copulaparam('t', 0.1, nu, 'type', 'spearman');
X2 = copularnd('t', rho_Pearson, nu, N);
figure(); scatterhist(X2(:,1),X2(:,2)); title('t, nu=1');

Consider two random variables X and Y both uniformly distributed in [0,1] and correlated with correlation rho.

I want to draw P realisations from their joint distribution in Matlab.

Is

A = copularnd('gaussian',rho,P);

the only way to do it?

解决方案

Copulas provide a very convenient way modelling the joint distribution. Just saying that you want X~U[0,1], Y~U[0,1] and corr(X,Y)=rho is not enough to define the relationship between these two random variables. Simply by substituting different copulas you can build different models (eventually choosing the one suiting your use-case best) that all satisfy this condition.

Apart from understanding the basics of what copulas are, you need to understand there are different types of correlations, such as linear (Pearson) and rank (Spearman / Kendall), and how they relate to each other. In particular, rank correlations preserve when a monotonic transformation is applied, unlike linear correlation. This is the key property allowing us to easily translate the desired linear correlation of uniform marginal distributions to the linear correlation of bivariate normal (or other type of distribution you use in the copula), which would be the input correlation to copularnd.

There is a very good answer on Cross Validated, describing exactly how to convert the correlation in your case. You should also read this Matlab guide on Using Rank Correlation Coefficients and copulas.

In a nutshell, to model desired linear correlation of marginals, you need to translate it into some rank correlation (using Spearman is convenient, since for uniform distribution it equals Pearason correlation), which would be the same for your normals (because it's a monotonic transformation). All you need is to convert that Spearman correlation for your normal to a linear Person correlation, which would be the input to copularnd.

Hopefully the code would make it easy to understand. Obviously you don't need those temporary variables, but it should make the logic clear:

rho_uni_Pearson = 0.7;
rho_uni_Spearman = rho_uni_Pearson;
rho_normal_Spearman = rho_uni_Spearman;
rho_normal_Pearson = 2*sin(rho_normal_Spearman*pi/6);

X = copularnd('gaussian', rho_normal_Pearson, 1e7);

And the resulting linear correlation is exactly what we wanted (since we've generated a very large sample):

corr(X(:,1), X(:,2));
ans =
    0.7000

Note that for a bivariate normal copula, the relationship between linear and rank correlations is easy enough, but it could be more complex if you use a different copula. Matlab has a function copulaparam that allows you to translate from rank to linear correlation. So instead of writing out explicit formula as above, we could just use:

rho_normal_Pearson = copulaparam('Gaussian', rho_normal_Spearman, 'type', 'spearman')

Now that we have learned the basics, let's go ahead and use a t copula with 5 degress of freedom instead of gaussian copula:

nu = 5; % number of degrees of freedom of t distribution
rho_t_Pearson = copulaparam('t', 0.7, nu, 'type', 'spearman');

X = copularnd('t', rho_t_Pearson, nu, 1e7);

Ensuring resulting linear correlation is what we wanted:

corr(X(:,1), X(:,2));
ans =
    0.6996

It's easy to observe that resulting distributions may be strikingly different, depending on which copula you choose, even through they all give you the same linear correlation. It is up to you, the researcher, to determine which model is best suited to your particular data and problem. For example:

N = 500;
rho_Pearson = copulaparam('gaussian', 0.1, 'type', 'spearman');
X1 = copularnd('gaussian', rho_Pearson, N);
figure(); scatterhist(X1(:,1),X1(:,2)); title('gaussian');

nu = 1; % number of degrees of freedom of t distribution
rho_Pearson = copulaparam('t', 0.1, nu, 'type', 'spearman');
X2 = copularnd('t', rho_Pearson, nu, N);
figure(); scatterhist(X2(:,1),X2(:,2)); title('t, nu=1');

这篇关于从Matlab中相关的边际统一中提取:copularnd的替代方案?的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆