A fast svd implementation suits your needs. The first place to start is to try scipy.sparse.linalg.svds. If it's not fast enough for your needs, then you'll need to explore other implementations. Some are mentioned here What fast algorithms exist for computing truncated SVD? but you might have to do some digging to turn up a python implementation that suits your needs. I've read that implicitly restarted lanczos bidiagonalization methods can compute SVD for much, much larger matrices than the one you have, such as the Netflix prize dataset, but I haven't found a widely used python implementation.
Alternatively, you could simply use your rank-deficient matrix and estimate a penalized regression instead. Options include ridge, lasso and elastic-net.