Matlab code for generating random symmetric positive definite matrix

Could anybody tell me How to generate random symmetric positive definite matrix using Matlab? Thank you very much for the help and suggestions.

-
Random with which distribution? –  Lord Soth Apr 11 at 5:46
@LordSoth Uniform distribution –  srijan Apr 11 at 5:47
I can give an algorithm that will generate a "random" symmetric, positive definite matrix, but the entries are by no means uniformly distributed, if they follow a standard distribution at all. –  Daryl Apr 11 at 6:08
The set of symmetric positive definite matrices is not compact, so such a thing as uniform distribution does not exist. –  Jyrki Lahtonen Apr 11 at 9:59
Whoever tells you to do that should then also specify the distribution. If not directly, then via a description of the random process that you are expected to study. I suspect that Wishart would be good one (see Johnny's answer). But really your task has not been fully specified, so your responsibility might be to go to your boss, and ask for more information - informing him/her about the danger of "garbage in/ garbage out" simulations to be done otherwise :-) –  Jyrki Lahtonen Apr 11 at 10:10

The algorithm I described in the comments is elaborated below. I will use $\tt{MATLAB}$ notation.

function A = generateSPDmatrix(n)
% Generate a dense n x n symmetric, positive definite matrix

A = rand(n,n); % generate a random n x n matrix

% construct a symmetric matrix using either
A = A+A'; OR
A = A*A';
% The first is significantly faster: O(n^2) compared to O(n^3)

% since A(i,j) < 1 by construction and a symmetric diagonally dominant matrix
%   is symmetric positive definite, which can be ensured by adding nI
A = A + n*eye(n);

end


Several changes are able to be used in the case of a sparse matrix.

function A = generatesparseSPDmatrix(n,density)
% Generate a sparse n x n symmetric, positive definite matrix with
%   approximately density*n*n non zeros

A = sprandsym(n,density); % generate a random n x n matrix

% since A(i,j) < 1 by construction and a symmetric diagonally dominant matrix
%   is symmetric positive definite, which can be ensured by adding nI
A = A + n*speye(n);

end


In fact, if the desired eigenvalues of the random are known and stored in rc, then the command

A = sprandsym(n,density,rc);


will construct the desired matrix. (Source: MATLAB sprandsym website)

-
Daryl Thank you very much daryl. It worked for me. Many thanks. –  srijan Apr 11 at 9:28
@Daryl, Thanks a lot! –  nullgeppetto Nov 20 at 9:39

A usual way in Bayesian statistics is to sample from a probability measure on real symmetric positive-definite matrices such as Wishart (or Inverse-Wishart).

I don't use Matlab but a quick check on Google gives this command (available in the Statistics toolbox):

W = wishrnd(Sigma,df)

where Sigma is some user-fixed positive definite matrix such as the identity and df are degrees of freedom.

-
Dear I need matlab code. I am not using any statistical software. –  srijan Apr 11 at 8:24
That is Matlab code! But you may or may not have Statistics toolbox. Try running 'wishrnd(eye(10),10)'. If it works, good for you. –  johnny Apr 11 at 8:28
Thank you very much. This code is working for me. I have one last question can we normalize the entries of the matrix to make them lie on some particular interval say $[0, 1]$ –  srijan Apr 11 at 9:05
The following is not computationally efficient but very simple. You could fill a matrix $\bf A$ with random values, computed for some desired distribution. Then you define a new matrix $\bf B = \bf{A} + \bf{A}^T$ in order to get a symmetric matrix. Then you use matlab to compute the eigenvalues of this matrix. Finally, you construct your matrix by
$$\bf{C} = \bf{B} + (\lambda_{min} + \delta)\bf{I}$$
where $\lambda_{min}$ is the smallest eigenvalue of $\bf{B}$ and $\delta$ is some small positive constant which defines the smallest eigenvalue of the your final matrix $\bf{C}$.
Matrix $\bf{B}$ is only symmetric, but not necessarily positive definite. –  Matt L. Apr 11 at 8:32
How to find $\lambda_{min}$ in case of random matrices? Eigenvalues will keep on varying every time I run the program. –  srijan Apr 11 at 8:43