Take the 2-minute tour ×
Cross Validated is a question and answer site for people interested in statistics, machine learning, data analysis, data mining, and data visualization. It's 100% free, no registration required.

My goal is to generate two variable which are correlated and one of them is heteroscedastic with regards to an grouping variable.

To create two variables with a desired correlation the common way to go is to use the cholesky decomposition. I additionally used the residuals of the orthogonal projection and standardized the variables up front to make sure the correlation stays at the desired value.

To introduce heteroscedasticity to one of the variables i tried the following:

  1. I did the cholesky decomposition for $u=(u_1,u_2), \ q[i] = c(1,1.01,1.02,1.03,...)$

    $Var(u) = \begin{pmatrix} 1 & \rho\cdot q[i] \\ \rho\cdot q[i] & q[i]^2 \end{pmatrix}$

    The reason why I choose this covariance matrix is that while the variance of $u_2$ increases with $q[i]$ the correlation stays at the desired value $\rho$.

    But this does obviously not introduce heteroscedasticity to $u_2$....

  2. A second attempt was to decomposed both variables $u_1,u_2$ with sample size n into two parts each of length $\frac{1}{2}\cdot n$: $\ u_1 = (u_{1.1},u_{1.2}), u_2 = (u_{2.1},u_{2.2})$. Now i did the cholesky decomposition for

    $Var(u) = \begin{pmatrix} 1 & 0 & \rho & \rho\sqrt{q[i]} \\ 0 & 1 & \rho & \rho\sqrt{q[i]} \\ \rho & \rho & 1 & 0 \\ \rho\sqrt{q[i]} & \rho\sqrt{q[i]} & 0 & q[i] \end{pmatrix}$

    After I generated the four variables I merged $u_{1.1},u_{1.2}$ back to one $n\times 1$ vector. Analogous for $u_2$. Testing the outcome yields that the heteroscedasticity is now where it should be as the variance of $u_{2.2}$ increases with $q[i]$. But the correlation vanishes as $q[i]$ increases because the variance of $u_2$ increases more rapidly then the covariance between $u_1$ and $u_2$

  3. A third attempt which i end up doing was to just multiply $u_2$ with a dummy-variable

    $z_2 =\begin{cases} 1\\ q[i] \end{cases}$

    This approach does work for small values of $q[i]$ (as my graphic shows) but as soon as $q[i]$ is large enough it dominates the nominator and denominator of the correlation coefficient such that the conditional variance (depending on the value of $z_2$) does no longer increase in $q[i]$.

Here is my code for case (3):

q     <- seq(1.01,10,0.1)
n         <- 100
rho   <- 0.5
sd_u1     <- numeric(0)
sd_u2.1   <- numeric(0)
sd_u2.2   <- numeric(0)
cor_u1_u2 <- numeric(0)


for(i in 1:length(q)){
u1      <- rnorm(n,0,1)
u1      <- ( u1 - mean(u1) )/sd(u1)
u2      <- rnorm(n,0,1)
z2      <- c(rep(1,0.5*n),rep(q[i],0.5*n))
u2      <- u2*z2
u2      <- as.vector( ( diag(n) - u1%*%solve(t(u1)%*%u1)%*%t(u1) ) %*% u2 )
u2      <- ( u2 - mean(u2) )/sd(u2)
z       <- cbind(u1, rho*u1+sqrt(1-rho^2)*u2) 
sd_u1[i]      <- sd(z[,1])
sd_u2.1[i]    <- sd(z[,2][1:(0.5*n)])
sd_u2.2[i]    <- sd(z[,2][(0.5*n+1):n])
cor_u1_u2[i]  <- cor(z[,1],z[,2])
}
par(mfrow=c(3,1))
plot(q,sd_u1, type="l")
plot(q,sd_u2.1, type="l", ylim=c(0,2))
lines(q,sd_u2.2,col="red")
plot(q,cor_u1_u2, type="l")

enter image description here

share|improve this question
    
Your code doesn't run as is. Please post working code. –  Max Nov 8 '12 at 3:35
    
now it should run. thx btw –  Druss2k Nov 8 '12 at 4:50
2  
It would help to explain your problem in English and mathematical terms, rather than wholly as an R program. There are several reasons for this, including (1) the program might not be doing what you intend, (2) it would help readers avoid having to work through code to understand the question, and (3) in the course of having to explain the circumstances clearly, you might arrive at a solution yourself. –  whuber Nov 8 '12 at 14:59
    
i edited a short description on top of the main article –  Druss2k Nov 8 '12 at 16:34
1  
Let's check that I understand, because many interpretations still seem possible. Do you wish to generate a sequence of vector-valued variables $(U_i,V_i)$ such that (1) Var$(V_i)=q_i$ (with $q_i$ predetermined) and (2) the correlation coefficient of the realized values $(u_i,v_i)$ equals $\rho$? (Please notice the distinction between the random variables and their realizations in this description.) And perhaps you want $V_i/q_i$ and $U_i$ all to be identically distributed? Maybe you only want the expectation of the correlation coefficient of $(u_i,v_i),i=1,\ldots,n$ to equal $\rho$? –  whuber Nov 8 '12 at 20:35

1 Answer 1

up vote 1 down vote accepted

Ruscio and Kaczetow (2008) proposed a method for inducing an exact correlation into two random variables with arbitrary distributions. Both variables can be sampled from existing empirical distributions, and the marginal distributions will reflect these distributions.

They have an accompanying R code which works (I used it a couple of times). As some minor typing errors in the original code had to be corrected: here's a cleaned version.

Ruscio, J., & Kaczetow, W. (2008). Simulating Multivariate Nonnormal Data Using an Iterative Algorithm. Multivariate Behavioral Research, 43(3), 355–381. doi:10.1080/00273170802285693

share|improve this answer
    
I hope that is an answer to the question (as it is not entirely clear to me as well ...) –  Felix S Nov 9 '12 at 13:46
    
sry i was at my moms :) yes this was very helpfull though i did it in a different way. as i said i was looking for a way to generate two variables with a certain correlation from which one is heteroscedastic. one covariance-matrix which does the job is the following: $\begin{pmatrix} 1 & 0 & \rho & rho[l]*\sqrt{q[i]} \\ 0 & 1 & \rho & \rho*\sqrt{2}*\sqrt{1+q[i]}-\rho \\ \rho & \rho & 1 & 0 \\ \rho*\sqrt{q[i]} & \rho*\sqrt{2}*\sqrt{1+q[i]}-\rho & 0 & q[i] \end{pmatrix}$ –  Druss2k Nov 11 '12 at 14:01
    
Before I was using the Cholesky-Decomposition to get the desired weights but I find it much easier to use the Singular Value Decomposition instead. Though this matrix does control the correlation in each group it tends to get larger in the group $u_{1.2},u_{2.2}$ for a certain $\rho$ if $q[i]$ increases. If $\rho$ reaches 0.5 the overall correlation does does decrease too if $q[i]$ increases. I know that my descriptions are rather cryptic and im very sorry about that. If someone is really interestend in my simulation study i would really like to explained it in more detial perhaps via chat :) –  Druss2k Nov 11 '12 at 14:06

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.