I have a python function that follows:
@njit(cache=True, fastmath=True, parallel=True)
def min_var_opt(covariance_matrix: np.ndarray) -> np.array:
n = covariance_matrix.shape[0]
# Compute the inverse of the covariance matrix
inverse_covariance = np.linalg.inv(covariance_matrix)
ones_vector = np.ones((n, 1))
denominator = np.dot(np.dot(ones_vector.T, inverse_covariance), ones_vector)
weights = np.dot(inverse_covariance, ones_vector) / denominator
return weights
covariance_matrix = np.array(
[
[77.25607201, 4.95549419, -2.1582784],
[4.95549419, 73.87388998, -3.76609601],
[-2.1582784, -3.76609601, 259.46734795],
]
)
weights = min_var_opt(covariance_matrix).T
print(f"compiled output:\t{weights}")
This produces the correct output: compiled output: [[0.41740475 0.4411879 0.14140735]]
Now, I would like to generalize this over a series of covariance matrices:
@guvectorize([(float64[:, :], float64[:])], "(n, n) -> (n)", nopython=True, cache=True)
def minimum_variance_optimization(covariance_matrix: np.ndarray, weights: np.array):
n = covariance_matrix.shape[0]
# Compute the inverse of the covariance matrix
inverse_covariance = np.linalg.inv(covariance_matrix)
ones_vector = np.ones((n, 1))
denominator = np.dot(np.dot(ones_vector.T, inverse_covariance), ones_vector)
weights = np.dot(inverse_covariance, ones_vector) / denominator
cov = np.array([covariance_matrix])
vectorized_weights = minimum_variance_optimization(cov)
print(f"vectorized output once:\t{vectorized_weights}")
This gives the identical output: vectorized output once: [[0.41740475 0.4411879 0.14140735]]
This is the expected output. The issue seems to be when i apply it with more than once inner matrix.
When I try:
cov = np.tile(covariance_matrix, (5, 1, 1))
vectorized_weights = minimum_variance_optimization(cov)
print(f"vectorized output repeated:\t{vectorized_weights}")
I get the following incorrect output:
vectorized output repeated:[[4.66464261e-310 0.00000000e+000 0.00000000e+000]
[4.66394411e-310 3.95252517e-323 0.00000000e+000]
[3.14015482e+179 2.36897524e+261 2.93013167e-002]
[7.10264216e+083 1.32448107e+007 3.12881791e+011]
[6.73475158e+107 4.17959458e+025 5.71862049e-310]]
According to this answer, I need to change the assignment above to be weights[:] = np.dot(inverse_covariance, ones_vector) / denominator to override the input array, but when I try this, numba fails to compile, and I don't understand how or why but it is writing a result with the original code above
I expect the output to be the same as the first vector repeated 5 times. I must be invoking the the gufunc incorrectly, but I cannot find any examples that match this use case. Can someone help me understand what is happening?
You need to use for loop for the first dimension of input array. Here is fixed implementation:
@nb.guvectorize([(nb.float64[:, :, :], nb.float64[:, :])], "(k, n, n) -> (k, n)", nopython=True, cache=True)
def minimum_variance_optimization(covariance_matrix: np.ndarray, weights: np.ndarray):
for k in range(covariance_matrix.shape[0]):
n = covariance_matrix.shape[1]
# Compute the inverse of the covariance matrix
inverse_covariance = np.linalg.inv(covariance_matrix[k, ])
ones_vector = np.ones((n, 1))
denominator = np.dot(np.dot(ones_vector.T, inverse_covariance), ones_vector)
weights[k] = (np.dot(inverse_covariance, ones_vector) / denominator)[:, 0]
Results:
cov = covariance_matrix[np.newaxis, ]
vectorized_weights = minimum_variance_optimization(cov)
print(f"vectorized output once:\t{vectorized_weights}")
vectorized output once: [[0.41740475 0.4411879 0.14140735]]
cov = np.tile(covariance_matrix, (5, 1, 1))
vectorized_weights = minimum_variance_optimization(cov)
print(f"vectorized output repeated:\t{vectorized_weights}")
vectorized output repeated: [[0.41740475 0.4411879 0.14140735]
[0.41740475 0.4411879 0.14140735]
[0.41740475 0.4411879 0.14140735]
[0.41740475 0.4411879 0.14140735]
[0.41740475 0.4411879 0.14140735]]