After finding a few unanswered requests for a solution of this problem in the web (including my own…) I’d like to share the final results of my work.

The problem:

Suppose you have two random variables, Z and T.

Z is N(0,1) distributed.

T is t(3) distributed.

Now you are supposed to produce four contour plots of the random variables’ joint pdf for the cases that the variables’ dependence structure is given by the

- Gaussian,
- Clayton,
- Frank- and
- Gumbel copula.

With the copula and the marginal distributions given the (bivariate) joint distribution of Z and T can be constructed. And this post is about doing exactly this in R and MatLab (and drawing the corresponding contour-plots).

### Implementation in R

The function `mvdc`

of the copula-package makes the solution in R quite easy:

library(copula)

layout(matrix(c(1,2,4,3), 2, 2, byrow = TRUE))

par(pty=”s”)#gaussian copula

gaussMVD<-mvdc(normalCopula(0), margins=c(“norm”,”t”), paramMargins=list(list(mean=0,sd=1), list(df=3)))

contour(gaussMVD, dmvdc, xlim=c(-3, 3), ylim=c(-3, 3), nlevels=10, xlab=”X”, ylab=”Y”, cex.axis = 1.5)#clayton copula

clayMVD<-mvdc(claytonCopula(0.9), margins=c(“norm”,”t”), paramMargins=list(list(mean=0,sd=1), list(df=3)))

contour(clayMVD, dmvdc, xlim = c(-3, 3), ylim = c(-3, 3), xlab=”X”, ylab=”Y”, cex.axis = 1.5)#gumbel copula

gumMVD<-mvdc(gumbelCopula(2), margins=c(“norm”,”t”), paramMargins=list(list(mean=0,sd=1), list(df=3)))

contour(gumMVD, dmvdc, xlim = c(-3, 3), ylim = c(-3, 3), xlab=”X”, ylab=”Y”, cex.axis = 1.5)#frank copula

frankMVD<-mvdc(frankCopula(8), margins=c(“norm”,”t”), paramMargins=list(list(mean=0,sd=1), list(df=3)))

contour(frankMVD, dmvdc, xlim = c(-3, 3), ylim = c(-3, 3), xlab=”X”, ylab=”Y”, nlevels=10, cex.axis = 1.5)

The result should look like this:

*God bless the incredible amount of packages available for R.*

### Implementation in MatLab

Solving the problem in MatLab is a little more tricky because there is no function like `mvdc`

available as is in R – at least as I know.

Of course – as in most cases – you could improve the code but as loops are very slow in MatLab I didn’t even try to use any. Anyway I’m an economist and business administrator and not a programmer. So I’m always happy to have managed converting some sophisticated problem into a working program (which no quantitative finance guy I asked has been able to help me with – just to have stated it…) .

Okay, here’s my code:

%with Gaussian copula

p = 0; % parameter of Gaussian copula

u=linspace(-3,3,1000); % 1000 equidistant numbers

[x,y]=meshgrid(u); % create the grid (the points) on which the functions are evaluated

u = normcdf(x,0,1); % compute standard normal cdf values for x

v = tcdf(y,3); % compute t3 cdf values for y

a = normpdf(x,0,1); % compute standard normal pdf values for x

b = tpdf(y,3); % compute t3 pdf values for y

f = copulapdf(‘Gaussian’,[u(:),v(:)],p); %compute the pdf of gaussian copula

mvd = reshape(f,1000,1000).*a.*b; % create multivariate density from the product

% of the copula pdf and the marginal pdfs

% reshape f to avoid matrix

% dimension conflictsfigure;

subplot(1,2,1);

contour(x,y,mvd);

title(‘Gaussian(0) Copula’);

axis square;

subplot(1,2,2);

mesh(x,y,mvd);

axis square;%———————————————————-

%with Clayton copulap = 0.9; % parameter of Clayton copula

u=linspace(-4,4,1000);

[x,y]=meshgrid(u);

u = normcdf(x,0,1);

v = tcdf(y,3);

a = normpdf(x,0,1);

b = tpdf(y,3);

f = copulapdf(‘Clayton’,[u(:),v(:)],p);

mvd = reshape(f,1000,1000).*a.*b;figure;

subplot(1,2,1);

contour(x,y,mvd);

title(‘Clayton(0.9) Copula’);

axis square;

subplot(1,2,2);

mesh(x,y,mvd);%——————————————————-

%with Frank copulap = 8; % parameter of Frank copula

u=linspace(-4,4,1000);

[x,y]=meshgrid(u);

u = normcdf(x,0,1);

v = tcdf(y,3);

a = normpdf(x,0,1);

b = tpdf(y,3);

f = copulapdf(‘Frank’,[u(:),v(:)],p);

mvd = reshape(f,1000,1000).*a.*b;figure;

subplot(1,2,1);

contour(x,y,mvd);

title(‘Frank(8) Copula’);

axis square;

subplot(1,2,2);

mesh(x,y,mvd);%%%%%%%%%%%%%%%%%%%%%%%%

%%with Gumbel copulap = 2; % parameter of Gumbel copula

u=linspace(-4,4,1000);

[x,y]=meshgrid(u);

u = normcdf(x,0,1);

v = tcdf(y,3);

a = normpdf(x,0,1);

b = tpdf(y,3);

f = copulapdf(‘Gumbel’,[u(:),v(:)],p);

mvd = reshape(f,1000,1000).*a.*b;figure;

subplot(1,2,1);

contour(x,y,mvd);

title(‘Gumbel(2) Copula’);

axis square;

subplot(1,2,2);

mesh(x,y,mvd);

Yuri Salazar FlroesSuppose I want the same but making no assumption on the marginals how would I do it? I tried assuming uniform margins but didn’t work, thanks in advance

adminHi Yuri.

You need to assume specific marginals, because otherwise one “variable” is missing to construct the bivariate pdf. Assuming uniform marginals is such an assumption.

What exactly didn’t work out for you? If you have problems with the program code, please note that the quotation marks may cause problems when copying them from the blog post and pasting them into the programs. You may need to replace them by your keyboard characters ” or ‘.

Best regards,

Christian

BenHi Christian

I have used this code for generating bivariate but i have a problem i have to generate pdf for mixing probabilities. Can you tell me how can I give input variable based on two different probability distribution.

Thanks in advance

Beena