[R] how to generate a random correlation matrix with restrictions

[Kingsford Jones]

Thanks Ravi -- with my limited linear algebra skills I was assuming
invertible symmetric was sufficient (rather than just necessary) for
positive definiteness.  So, the open challenge is to generate a pd
matrix of dimension 100 with r_ij = 1 if i=j else -1< r_ij< 1, with
about 10% of the elements >.9.

I tried for a bit without luck, but I can offer another way to
generate a pd matrix:

ev <- runif(100) # choose some positive eigenvalues
U <- svd(matrix(runif(10000), nc=100))$u  # an orthogonal matrix
pdM <- U %*% diag(ev) %*% t(U)


Kingsford



On Fri, Aug 28, 2009 at 9:41 PM, Ravi Varadhan<rvaradhan at jhmi.edu> wrote:
> Hi Kingsford,
>
> There is more structure to a correlation matrix than that meets the eye!  Your method will produce a matrix R that looks "like" a correlation matrix, but beware - it is an impostor!
>
> You can obtain a valid correlation matrix, Q, from the impostor R by using the `nearPD' function in the "Matrix" package, which finds the positive definite matrix Q that is "nearest" to R.  However, note that when R is far from a positive-definite matrix, this step may give a Q that does not have the desired property.
>
> require(Matrix)
>
> R <- matrix(runif(16), ncol=4)
>
> R <- (R * lower.tri(R)) + t(R * lower.tri(R))
>
> diag(R) <- 1
>
> eigen(R)$val
>
> Q <- nearPD(R, posd.tol=1.e-04)$mat
>
> eigen(Q)$val
>
> max(abs(Q - R))  # maximum discrepancy between R and Q
>
> Another easy way to produce a valid correlation matrix is:
>
> R <- matrix(runif(36), ncol=6)
>
> RtR <- R %*% t(R)
>
> Q <- cov2cor(RtR)
>
> But this does not have the property that the correlations are uniformly distributed.
>
> Hope this helps,
> Ravi.
>
> ____________________________________________________________________
>
> Ravi Varadhan, Ph.D.
> Assistant Professor,
> Division of Geriatric Medicine and Gerontology
> School of Medicine
> Johns Hopkins University
>
> Ph. (410) 502-2619
> email: rvaradhan at jhmi.edu
>
>
> ----- Original Message -----
> From: Kingsford Jones <kingsfordjones at gmail.com>
> Date: Friday, August 28, 2009 10:12 pm
> Subject: Re: [R] how to generate a random correlation matrix with       restrictions
> To: Ning Ma <pningma at gmail.com>
> Cc: r-help at r-project.org
>
>
>> Ahh -- Mark Leeds just pointed out off-list that it was supposed to be
>> a correlation matrix.
>>
>> Perhaps
>>
>> R <- matrix(runif(10000), ncol=100)
>> R <- (R * lower.tri(R)) + t(R * lower.tri(R))
>> diag(R) <- 1
>>
>> These are of course uniformly distributed positive correlations.
>>
>> Kingsford
>>
>>
>> On Fri, Aug 28, 2009 at 7:36 PM, Kingsford
>> Jones<kingsfordjones at gmail.com> wrote:
>> >  R <- matrix(runif(10000), ncol=100)
>> >
>> > hth,
>> >
>> > Kingsford Jones
>> >
>> > On Fri, Aug 28, 2009 at 7:26 PM, Ning Ma<pningma at gmail.com> wrote:
>> >> Hi,
>> >>
>> >> How can I generate a random 100x100 correlation matrix, R={r_ij},
>> >> where about 10% of r_ij are greater than 0.9
>> >>
>> >> Thanks in advance.
>> >>
>> >> ______________________________________________
>> >> R-help at r-project.org mailing list
>> >>
>> >> PLEASE do read the posting guide
>> >> and provide commented, minimal, self-contained, reproducible code.
>> >>
>> >
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>>
>> PLEASE do read the posting guide
>> and provide commented, minimal, self-contained, reproducible code.
>