Generating a sparse matrix from a piecewise function

General Usage Performance
, ,
1

I’m trying to generate a random sparse matrix R where each element in R can take on the following values:

r_{ij}= \left\{ \begin{array}{ll} 1 & \text{with probability } \frac{1}{6}\\ 0 & \text{with probability } \frac{2}{3} \\ -1 & \text{with probability } \frac{1}{6} \\ \end{array} \right.

What would be the best way to do this in Julia? My first approach looks like this:

using SparseArrays

function genR1(n,m)
    R = spzeros(n,m)

    for i in eachindex(R)
        r = rand()
        if r > 5/6
            R[i] = 1
        elseif r > 2/3
            R[i] = -1
        end
    end
    
    return R
end

however, I found this to be slower and allocate more than not using sparse matrices at all.

2

Don’t use a sparse matrix for this. Sparse matrices are only for sparsity of 95% or greater.

2 Likes
3

No matter what you do, mutating a sparse array element-by-element is very slow (because the compressed column format needs to be re-packed each time you insert a new nonzero).

Here, your matrix isn’t very sparse, so I wouldn’t use a sparse matrix at all in this case. You could just do

R = rand([1,-1,0,0,0,0], n, m)

If you want you can convert this to a sparse matrix with sparse(R), but it’s generally not worth using a sparse-matrix data structure with a matrix that is only 2/3 sparse.

However, it is possible to generate such a sparse matrix directly, by calling sprand with a custom random-number generator. (This would be worth it if the probability p of a nonzero entry were much smaller than 1/3.) In particular, you could do:

p = 1/3 # probability of nonzero (±1) entry
R = sprand(n, m, p, N -> rand((-1,1), N))
6 Likes