www.gusucode.com > stats 源码程序 matlab案例代码 > stats/SimulatingDependentMultivariateDataUsingCopulasExample.m

    %% Simulating Dependent Multivariate Data Using Copulas
%% Load sample data and plot the histograms

% Copyright 2015 The MathWorks, Inc.

load stockreturns
nobs = size(stocks,1);

subplot(2,1,1)
histogram(stocks(:,1),10,'FaceColor',[.8 .8 1])
xlim([-3.5 3.5])
xlabel('X1')
ylabel('Frequency')

subplot(2,1,2)
histogram(stocks(:,2),10,'FaceColor',[.8 .8 1])
xlim([-3.5 3.5])
xlabel('X2')
ylabel('Frequency')
%% Compare the empirical cdf to a kernel smoothed cdf estimate for the first variable
[Fi,xi] = ecdf(stocks(:,1));

figure()
stairs(xi,Fi,'b','LineWidth',2)
hold on

Fi_sm = ksdensity(stocks(:,1),xi,'function','cdf','width',.15);

plot(xi,Fi_sm,'r-','LineWidth',1.5)
xlabel('X1')
ylabel('Cumulative Probability')
legend('Empirical','Smoothed','Location','NW')
grid on
%% Compute the rank correlation of the data
nu = 5;
tau = corr(stocks(:,1),stocks(:,2),'type','kendall')
%% Find the corresponding linear correlation parameter for the t copula
rho = copulaparam('t', tau, nu, 'type','kendall')
%% Generate random values from the t copula
n = 1000;
U = copularnd('t',[1 rho; rho 1],nu,n);
%% Make a kernel estimate of distribution and evaluate the inverse cdf at the copula points 
X1 = ksdensity(stocks(:,1),U(:,1),...
               'function','icdf','width',.15);
X2 = ksdensity(stocks(:,2),U(:,2),...
               'function','icdf','width',.15);
%% compute the inverse cdf over a grid of values in the interval (0,1) and use interpolation to evaluate it at the copula points
p = linspace(0.00001,0.99999,1000);
G1 = ksdensity(stocks(:,1),p,'function','icdf','width',0.15);
X1 = interp1(p,G1,U(:,1),'spline');
G2 = ksdensity(stocks(:,2),p,'function','icdf','width',0.15);
X2 = interp1(p,G2,U(:,2),'spline');

scatterhist(X1,X2,'Direction','out')