matlab蒙特卡罗变量变换的分布，使⽤copula仿真相关随机变量仿真输⼊之间的相关性
蒙特卡罗模拟的设计决策之⼀是为随机输⼊选择概率分布。为每个变量选择⼀种分布往往很简单，但确定输⼊之间应该存在什么样的相关性却可能不那么简单。理想情况下，仿真的输⼊数据应反映要建模的实际数量之间已知的相关性。然⽽，判断仿真中的任何相关性时可以依据的信息可能很少或根本没有，在这种情况下，最好的做法是尝试不同的可能性，以确定模型的敏感性。
但是，当输⼊数据的分布不是标准多元分布时，很难真正⽣成具有相关性的随机输⼊。⽽且，有些标准多元分布只能对⾮常有限的⼏种相关性进⾏建模。我们可以始终将输⼊视为各⾃独⽴，这很简便，但并不总是合理的，有可能导致错误的结论。
例如，在有关⾦融风险的蒙特卡罗模拟中，可能有⼀些表⽰各种保险损失来源的随机输⼊。这些输⼊可以建模为对数正态随机变量。我们会想了解两两输⼊间的相关性对仿真结果有何影响。确实，从实际数据中可能已经知道，相同的随机条件可对两个损失来源都产⽣影响，如果在仿真中忽略这⼀点，可能会导致错误的结论。
独⽴的对数正态随机变量的仿真⾮常繁琐。最简单的⽅法是使⽤ lognrnd 函数。此处，我们将使⽤ mvnrnd 函数⽣成 n 对独⽴的正态随机变量，然后计算它们的幂。请注意，这⾥使⽤的协⽅差矩阵是对⾓矩阵，即 Z 的各列之间是独⽴的。
n = 1000;
sigma = .5;
SigmaInd = sigma.^2 .* [1 0; 0 1]
SigmaInd =
0.2500 0
0 0.2500
ZInd = mvnrnd([0 0], SigmaInd, n);
XInd = exp(ZInd);
plot(XInd(:,1),XInd(:,2),'.');
axis equal;
axis([0 5 0 5]);
xlabel('X1');
ylabel('X2');
⽣成具有相关性的⼆元对数正态随机变量也很简单，可使⽤⾮零项不在对⾓线上的协⽅差矩阵。
rho = .7;
SigmaDep = sigma.^2 .* [1 rho; rho 1]
SigmaDep =
0.2500 0.1750
0.1750 0.2500
ZDep = mvnrnd([0 0], SigmaDep, n);
XDep = exp(ZDep);
第⼆个散点图说明这两种⼆元分布之间的差异。
plot(XDep(:,1),XDep(:,2),'.');
axis equal;
axis([0 5 0 5]);
xlabel('X1');
ylabel('X2');
很明显，第⼆个数据集中较⼤的 X1 值更倾向于与较⼤的 X2 值关联，⼩值之间也是如此。这种相关性由基础⼆元正态分布的相关参数 rho 决定。从仿真中得出的结论可能很⼤程度上取决于⽣成的 X1 和 X2 是否存在相关性。
本⽰例中的⼆元对数正态分布是⼀个简单的解决⽅法，它当然可以泛化应⽤于更⾼维度以及边缘分布为其他对数正态分布的情况。也有⼀些多元分布，例如，多元 t 分布和 Dirichlet 分布，它们分别⽤于仿真相关的 t 随机变量和 beta 随机变量。但是，简单的多元分布并不多，⽽且仅适⽤于边缘分布都为同⼀族分布(甚⾄是完全相同的分布)的情形。在很多情况下，这可能会成为限制。
构造相关⼆元分布的更通⽤⽅法
尽管上⾯创建⼆元对数正态分布的构造很简单，但它说明了⼀种更普遍适⽤的⽅法。⾸先，我们从⼆元正态分布中⽣成值对。两两变量间存在统计相关性，⽽且每个变量都具有正态边缘分布。接下来，分别对每个变量进⾏变换(指数函数)，将边缘分布变成对数正态分布。变换后的变量仍具有统计相关性。
如果能够到⼀种合适的变换，就可以将此⽅法泛化，从⽽为其他边缘分布⽣成相关⼆元随机向量。实际上，确实存在构造这种变换的⼀般⽅法，尽管不像求幂那样简单。
根据定义，将正态 CDF(这⾥⽤ PHI 表⽰)应⽤于标准正态随机变量将⽣成在区间 [0,1] 内均匀分布的随机变量。为说明这⼀点，假设 Z 具有标准正态分布，则 U = PHI(Z) 时，CDF 为
Pr{U <= u0} = Pr{PHI(Z) <= u0} = Pr{Z <= PHI^(-1)(u0)} = u0,
⽽这就是 U(0,1) 随机变量的 CDF。通过为仿真的正态值和仿真的变换值绘制的直⽅图证实了这⼀点。
n = 1000;
z = normrnd(0,1,n,1);
hist(z,-3.75:.5:3.75);
xlim([-4 4]);
title('1000 Simulated N(0,1) Random Values');
xlabel('Z');
ylabel('Frequency');
u = normcdf(z);
hist(u,.05:.1:.95);
title('1000 Simulated N(0,1) Values Transformed to U(0,1)');
matlab生成随机数
xlabel('U');
ylabel('Frequency');
现在，借⽤⼀元随机数⽣成理论，将任意分布 F 的逆 CDF 应⽤于 U(0,1) 随机变量，可得到分布正好为 F 的随机变量。这就是逆变换法。其证明过程与上述正向情形的证明过程基本相反。另⼀个直⽅图说明了变换为 gamma 分布的结果。
x = gaminv(u,2,1);
hist(x,.25:.5:9.75);
title('1000 Simulated N(0,1) Values Transformed to Gamma(2,1)');
xlabel('X');
ylabel('Frequency');
这个两步变换⽅法可以应⽤于标准⼆元正态分布的每个变量，从⽽⽣成具有任意边缘分布的相关随机变量。由于变换分别作⽤于每个成分，因此⽣成的两个随机变量甚⾄不需要具有相同的边缘分布。变换定义为
Z = [Z1 Z2] ~ N([0 0],[1 rho; rho 1])
U = [PHI(Z1) PHI(Z2)]
X = [G1(U1) G2(U2)]
其中 G1 和 G2 是两个可能不相同的分布的逆 CDF。例如，我们可以从具有 t(5) 和 Gamma(2,1) 边缘分布的⼆元分布⽣成随机向量。
n = 1000;
rho = .7;
Z = mvnrnd([0 0], [1 rho; rho 1], n);
U = normcdf(Z);
X = [gaminv(U(:,1),2,1) tinv(U(:,2),5)];
在下图中，直⽅图显⽰在散点图的左侧和下⽅，以显⽰边缘分布和相关性。
[n1,ctr1] = hist(X(:,1),20);
[n2,ctr2] = hist(X(:,2),20);
subplot(2,2,2);
plot(X(:,1),X(:,2),'.');
axis([0 12 -8 8]);
h1 = gca;
title('1000 Simulated Dependent t and Gamma Values');
xlabel('X1 ~ Gamma(2,1)');
ylabel('X2 ~ t(5)');
subplot(2,2,4);
bar(ctr1,-n1,1);
axis([0 12 -max(n1)*1.1 0]);
axis('off');
h2 = gca;
subplot(2,2,1);
barh(ctr2,-n2,1);
axis([-max(n2)*1.1 0 -8 8]);
axis('off');
h3 = gca;
h1.Position = [0.35 0.35 0.55 0.55];
h2.Position = [.35 .1 .55 .15];
h3.Position = [.1 .35 .15 .55];
colormap([.8 .8 1]);
秩相关系数
此构造中 X1 和 X2 之间的相关性由基础⼆元正态分布的相关参数 rho 决定。但是，X1 和 X2 的线性相关系数并不完全等于 rho。例如，在原来的对数正态分布中，该相关系数的闭式解为：
cor(X1,X2) = (exp(rho.*sigma.^2) - 1) ./ (exp(sigma.^2) - 1)
它严格⼩于 rho，除⾮ rho 刚好为 1。但更普遍的情形是，例如在上⾯的 gamma/t 分布中，X1 和 X2 之间的线性相关性很难或⽆法⽤rho 表⽰，但可以通过仿真来说明存在同样的效应。
这是因为线性相关系数表达的是随机变量之间的线性相关性，当对这些随机变量应⽤⾮线性变换后，将不会保留线性相关性。在这种情况下，秩相关系数(如 Kendall 的 tau 或 Spearman 的 rho)更合适。
⼤体上来讲，这些秩相关衡量的是⼀个随机变量的⼤值或⼩值与另⼀个随机变量的⼤值或⼩值之间的相关程度。然⽽，与线性相关系数不同的是，秩相关系数只从秩的⾓度衡量相关性。因此，秩相关在任何单调变换下都会保留。前⾯所述的变换⽅法也会保留秩相关。因此，知道⼆元正态分布 Z 的秩相关就可以准确地确定最后变换的随机变量 X 的秩相关。尽管仍然需要使⽤ rho 来参数化基础的⼆元正态分布，但Kendall tau 或 Spearman rho 在描述随机变量之间的相关性时更有⽤，因为它们不会随着选择的边缘分布不同⽽改变。

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。