Covariance¶
Note
This page describes the MaNGA DAP
Covariance
object in general
terms. For its specific use in the DAP for handling spatial
covariance in the MaNGA datacubes, see Spatial Covariance.
The DAP Covariance
object is a
general utility for constructing, visualizing, and storing two- and
three-dimensional covariance matrices.
Construction¶
Beyond the nominal instantiation method, there are numerous methods
used to construct a Covariance
object.
The simplest approaches are if you have a variance vector or a pre-calculated covariance matrix as a dense array and you want to construct a
Covariance
object for further calculations:import numpy from mangadap.util.covariance import Covariance # Construct the Covariance object just from a variance vector var = numpy.ones(3, dtype=float) covar = Covariance.from_variance(var) # Should be the same as the identity matrix. assert numpy.array_equal(covar.toarray(), numpy.identity(3)) # Construct the Covariance object from a pre-calculated array c = numpy.diag(numpy.full(10-2, 0.2, dtype=float), k=-2) \ + numpy.diag(numpy.full(10-1, 0.5, dtype=float), k=-1) \ + numpy.diag(numpy.full(10, 1.0, dtype=float), k=0) \ + numpy.diag(numpy.full(10-1, 0.5, dtype=float), k=1) \ + numpy.diag(numpy.full(10-2, 0.2, dtype=float), k=2) covar = Covariance.from_array(c) # Should be the same as the identity matrix. assert numpy.array_equal(covar.toarray(), c)Note the use of the
toarray()
to access the array; see Accessing the covariance data.You can construct a covariance matrix based on samples from a distribution:
import numpy from mangadap.util.covariance import Covariance # Build a bogus covariance matrix m = numpy.zeros(10, dtype=float) c = numpy.diag(numpy.full(10-2, 0.2, dtype=float), k=-2) \ + numpy.diag(numpy.full(10-1, 0.5, dtype=float), k=-1) \ + numpy.diag(numpy.full(10, 1.0, dtype=float), k=0) \ + numpy.diag(numpy.full(10-1, 0.5, dtype=float), k=1) \ + numpy.diag(numpy.full(10-2, 0.2, dtype=float), k=2) # Draw samples s = numpy.random.multivariate_normal(m, c, size=100000) # Instantiate covar = Covariance.from_samples(s.T, cov_tol=0.1) # Check the values are very nearly the same as the input assert numpy.all(numpy.absolute(c - covar.toarray()) < 0.02)Here, we’ve drawn samples from a known multivariate normal distribution with a given covariance between its 10 axes and checked the reconstruction of the covariance matrix from those samples using
from_samples()
. This construction method is a simple wrapper for numpy.cov andfrom_array()
.You can construct the covariance matrix that results from a matrix multiplication:
import numpy from mangadap.util.covariance import Covariance # Build a bogus covariance matrix c = numpy.diag(numpy.full(10-2, 0.2, dtype=float), k=-2) \ + numpy.diag(numpy.full(10-1, 0.5, dtype=float), k=-1) \ + numpy.diag(numpy.full(10, 1.0, dtype=float), k=0) \ + numpy.diag(numpy.full(10-1, 0.5, dtype=float), k=1) \ + numpy.diag(numpy.full(10-2, 0.2, dtype=float), k=2) # Vector and matrix to multiply x = numpy.ones(10, dtype=float) t = numpy.zeros((3,10), dtype=float) t[0,0] = 1.0 t[1,2] = 1.0 t[2,4] = 1.0 # Apply the matrix multiplication (not used; just for illustration) y = numpy.dot(t, x) # Expected covariance from the result _c = numpy.diag(numpy.full(3-1, 0.2, dtype=float), k=-1) \ + numpy.diag(numpy.full(3, 1.0, dtype=float), k=0) \ + numpy.diag(numpy.full(3-1, 0.2, dtype=float), k=1) covar = Covariance.from_matrix_multiplication(t, c) assert numpy.array_equal(covar.toarray(), _c)Finally, you can construct the covariance matrix from a previous instance that was saved to a fits file using the Fits file I/O methods.
Accessing the covariance data¶
The Covariance
object is primarily
a storage and IO utility. Internally, the object only keeps the upper
triangle of the matrix, which means that use of the cov
attribute is not recommended unless you know what you’re doing.
Also note that the object can be either “2D” or “3D”. When the object
is 3D (covar.dim == 3
), this is just a convenience for storing
multiple covariance matrices that have the same shape.
There are two ways to access the full covariance matrix, the
full()
and
toarray()
methods
depending on whether you want a sparse or dense matrix, respectively.
Also note that the use of the
full()
method requires you
to specify a single channel for 3D objects. The output of these two
methods can be used as you would use any scipy.sparse.csr_matrix
or numpy.ndarray object.
To show the covariance matrix, you can use its
show()
method to quickly
produce a plot, which is a simple wrapper for the
toarray()
method and
pyplot.imshow.
Toggling between covariance and correlation matrices¶
The Covariance
object allows you to
toggle between the full covariance matrix, \({\mathbf C}\) and a
correlation matrix, \({\mathbf \rho}\), where:
where \({\mathbf V}\) is the variance vector (the diagonal
elements of \({\mathbf C}\)). To convert a
Covariance
object to a correlation
matrix (or ensure that it already is one), use
to_correlation()
. To revert
back to a covariance matrix, use
revert_correlation()
.
Fits file I/O methods¶
Covariance
objects can be saved as
a binary table in a fits file using the
write()
method. To reload
the covariance matrix, use the
from_fits()
instantiation
method:
import numpy
from mangadap.util.covariance import Covariance
ofile = 'test_covar_io.fits'
# Build a bogus covariance matrix
m = numpy.zeros(10, dtype=float)
c = numpy.diag(numpy.full(10-2, 0.2, dtype=float), k=-2) \
+ numpy.diag(numpy.full(10-1, 0.5, dtype=float), k=-1) \
+ numpy.diag(numpy.full(10, 1.0, dtype=float), k=0) \
+ numpy.diag(numpy.full(10-1, 0.5, dtype=float), k=1) \
+ numpy.diag(numpy.full(10-2, 0.2, dtype=float), k=2)
# Draw samples
s = numpy.random.multivariate_normal(m, c, size=100000)
# Instantiate
covar = Covariance.from_samples(s.T, cov_tol=0.1)
# Write
covar.write(ofile)
# Read
_covar = Covariance.from_fits(ofile)
# Should be the same
assert numpy.allclose(covar.toarray(), _covar.toarray())
The details of how the covariance data are stored are described by
the write()
method.