<p dir="auto">Encountering a MemoryError when creating a large sparse random matrix.<br>
The matrix is created with scipy.sparse.rand, and even if the density is low enough that only a few elements in the whole matrix should be nonzero (or even none at all), it still tries to allocate a large integer array.</p>
<p dir="auto">I'm not sure if this really is a bug, or I'm just using the function wrong. I know it's possible to create an empty sparse matrix and fill it with random values manually. But I thought it a bit strange that it would try to allocate that much.</p>
<h4 dir="auto">Reproducing code example:</h4>

<div class="snippet-clipboard-content notranslate position-relative overflow-auto" data-snippet-clipboard-copy-content="from scipy.sparse import rand
rand(10**5, 10**5, density=0) "><pre class="notranslate"><code class="notranslate">from scipy.sparse import rand
rand(10**5, 10**5, density=0) 
</code></pre></div>
<h4 dir="auto">Error message:</h4>

<div class="snippet-clipboard-content notranslate position-relative overflow-auto" data-snippet-clipboard-copy-content="---------------------------------------------------------------------------
MemoryError                               Traceback (most recent call last)
&lt;ipython-input-68-4d3ba46fe0a7&gt; in &lt;module&gt;
----&gt; 1 rand(10**5, 10**5, density=0)

~/anaconda3/lib/python3.7/site-packages/scipy/sparse/construct.py in rand(m, n, density, format, dtype, random_state)
    840 
    841     &quot;&quot;&quot;
--&gt; 842     return random(m, n, density, format, dtype, random_state)

~/anaconda3/lib/python3.7/site-packages/scipy/sparse/construct.py in random(m, n, density, format, dtype, random_state, data_rvs)
    786             data_rvs = random_state.rand
    787 
--&gt; 788     ind = random_state.choice(mn, size=k, replace=False)
    789 
    790     j = np.floor(ind * 1. / m).astype(tp, copy=False)

mtrand.pyx in numpy.random.mtrand.RandomState.choice()

mtrand.pyx in numpy.random.mtrand.RandomState.permutation()

MemoryError: Unable to allocate 74.5 GiB for an array with shape (10000000000,) and data type int64
"><pre class="notranslate"><code class="notranslate">---------------------------------------------------------------------------
MemoryError                               Traceback (most recent call last)
&lt;ipython-input-68-4d3ba46fe0a7&gt; in &lt;module&gt;
----&gt; 1 rand(10**5, 10**5, density=0)

~/anaconda3/lib/python3.7/site-packages/scipy/sparse/construct.py in rand(m, n, density, format, dtype, random_state)
    840 
    841     """
--&gt; 842     return random(m, n, density, format, dtype, random_state)

~/anaconda3/lib/python3.7/site-packages/scipy/sparse/construct.py in random(m, n, density, format, dtype, random_state, data_rvs)
    786             data_rvs = random_state.rand
    787 
--&gt; 788     ind = random_state.choice(mn, size=k, replace=False)
    789 
    790     j = np.floor(ind * 1. / m).astype(tp, copy=False)

mtrand.pyx in numpy.random.mtrand.RandomState.choice()

mtrand.pyx in numpy.random.mtrand.RandomState.permutation()

MemoryError: Unable to allocate 74.5 GiB for an array with shape (10000000000,) and data type int64

</code></pre></div>
<h4 dir="auto">Scipy/Numpy/Python version information:</h4>
<div class="snippet-clipboard-content notranslate position-relative overflow-auto" data-snippet-clipboard-copy-content="1.4.1 1.18.1 sys.version_info(major=3, minor=7, micro=6, releaselevel='final', serial=0)"><pre class="notranslate"><code class="notranslate">1.4.1 1.18.1 sys.version_info(major=3, minor=7, micro=6, releaselevel='final', serial=0)
</code></pre></div>
<p dir="auto"><a class="issue-link js-issue-link" data-error-text="Failed to load title" data-id="341323117" data-permission-text="Title is private" data-url="https://github.com/scipy/scipy/issues/9036" data-hovercard-type="issue" data-hovercard-url="/scipy/scipy/issues/9036/hovercard" href="https://github.com/scipy/scipy/issues/9036">gh-9036</a> removed a special case in the <code class="notranslate">scipy.sparse.random</code> code, as the numpy implementation of <code class="notranslate">random.choice(..., replace=False)</code> had improved significantly.</p>
<p dir="auto">I think we may need to reintroduce the special case, but only for more extreme cases. A quick test shows that <code class="notranslate">np.random.choice(N, k, replace=False)</code> scales linearly in <code class="notranslate">N</code>:</p>
<div class="snippet-clipboard-content notranslate position-relative overflow-auto" data-snippet-clipboard-copy-content="In [2]: np.__version__
Out[2]: '1.15.4'

In [3]: %timeit np.random.choice(100, 10, replace=False)
17 µs ± 62.3 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

In [4]: %timeit np.random.choice(1000, 10, replace=False)
28.9 µs ± 204 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

In [5]: %timeit np.random.choice(10000, 10, replace=False)
152 µs ± 1.17 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

In [6]: %timeit np.random.choice(100000, 10, replace=False)
1.45 ms ± 6.24 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

In [7]: %timeit np.random.choice(1000000, 10, replace=False)
14.6 ms ± 101 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)"><pre class="notranslate"><code class="notranslate">In [2]: np.__version__
Out[2]: '1.15.4'

In [3]: %timeit np.random.choice(100, 10, replace=False)
17 µs ± 62.3 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

In [4]: %timeit np.random.choice(1000, 10, replace=False)
28.9 µs ± 204 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

In [5]: %timeit np.random.choice(10000, 10, replace=False)
152 µs ± 1.17 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

In [6]: %timeit np.random.choice(100000, 10, replace=False)
1.45 ms ± 6.24 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

In [7]: %timeit np.random.choice(1000000, 10, replace=False)
14.6 ms ± 101 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
</code></pre></div>
<p dir="auto">The result I've observed is that <code class="notranslate">scipy.sparse.random</code> gets unbearably slow. For example, when running the example code in <a class="issue-link js-issue-link" data-error-text="Failed to load title" data-id="400566780" data-permission-text="Title is private" data-url="https://github.com/scipy/scipy/issues/9697" data-hovercard-type="issue" data-hovercard-url="/scipy/scipy/issues/9697/hovercard" href="https://github.com/scipy/scipy/issues/9697">gh-9697</a> on scipy master, constructing the matrix takes over 2 minutes and completely locks my machine for the duration of the <code class="notranslate">np.random.choice</code> call.</p>
<p dir="auto">Perhaps the new logic should be:</p>
<div class="highlight highlight-source-python notranslate position-relative overflow-auto" dir="auto" data-snippet-clipboard-copy-content="    if mn &lt; 1e8 or k &gt; 1e4:
        ind = random_state.choice(mn, size=k, replace=False)
    else:
        ind = np.empty(k, dtype=tp)
        selected = set()
        for i in xrange(k):
            j = random_state.randint(mn)
            while j in selected:
                j = random_state.randint(mn)
            selected.add(j)
            ind[i] = j"><pre class="notranslate">    <span class="pl-k">if</span> <span class="pl-s1">mn</span> <span class="pl-c1">&lt;</span> <span class="pl-c1">1e8</span> <span class="pl-c1">or</span> <span class="pl-s1">k</span> <span class="pl-c1">&gt;</span> <span class="pl-c1">1e4</span>:
        <span class="pl-s1">ind</span> <span class="pl-c1">=</span> <span class="pl-s1">random_state</span>.<span class="pl-en">choice</span>(<span class="pl-s1">mn</span>, <span class="pl-s1">size</span><span class="pl-c1">=</span><span class="pl-s1">k</span>, <span class="pl-s1">replace</span><span class="pl-c1">=</span><span class="pl-c1">False</span>)
    <span class="pl-k">else</span>:
        <span class="pl-s1">ind</span> <span class="pl-c1">=</span> <span class="pl-s1">np</span>.<span class="pl-en">empty</span>(<span class="pl-s1">k</span>, <span class="pl-s1">dtype</span><span class="pl-c1">=</span><span class="pl-s1">tp</span>)
        <span class="pl-s1">selected</span> <span class="pl-c1">=</span> <span class="pl-en">set</span>()
        <span class="pl-k">for</span> <span class="pl-s1">i</span> <span class="pl-c1">in</span> <span class="pl-en">xrange</span>(<span class="pl-s1">k</span>):
            <span class="pl-s1">j</span> <span class="pl-c1">=</span> <span class="pl-s1">random_state</span>.<span class="pl-en">randint</span>(<span class="pl-s1">mn</span>)
            <span class="pl-k">while</span> <span class="pl-s1">j</span> <span class="pl-c1">in</span> <span class="pl-s1">selected</span>:
                <span class="pl-s1">j</span> <span class="pl-c1">=</span> <span class="pl-s1">random_state</span>.<span class="pl-en">randint</span>(<span class="pl-s1">mn</span>)
            <span class="pl-s1">selected</span>.<span class="pl-en">add</span>(<span class="pl-s1">j</span>)
            <span class="pl-s1">ind</span>[<span class="pl-s1">i</span>] <span class="pl-c1">=</span> <span class="pl-s1">j</span></pre></div>