-
-
Notifications
You must be signed in to change notification settings - Fork 18.5k
PERF: cythonize kendall correlation #39132
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Changes from 7 commits
72427e6
7ded648
ef1dc87
b2384bf
26f4771
e6fca9b
85a33c8
012b5a1
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -393,6 +393,99 @@ def nancorr_spearman(ndarray[float64_t, ndim=2] mat, Py_ssize_t minp=1) -> ndarr | |
return result | ||
|
||
|
||
# ---------------------------------------------------------------------- | ||
# Kendall correlation | ||
# Wikipedia article: https://en.wikipedia.org/wiki/Kendall_rank_correlation_coefficient | ||
|
||
@cython.boundscheck(False) | ||
@cython.wraparound(False) | ||
def nancorr_kendall(ndarray[float64_t, ndim=2] mat, Py_ssize_t minp=1) -> ndarray: | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. we usually dont mix-and-match cython-vs-python style annotations. i know this is tough bc the python-style won't allow There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Right now, the function signature for this function is the same as the one for nancorr_spearman. However, I notice that other functions have left out the There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. right but we are currently not checking the .pyx right so this is moot for now? There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. yah its not a big deal |
||
""" | ||
Perform kendall correlation on a 2d array | ||
|
||
Parameters | ||
---------- | ||
mat : Array to compute kendall correlation on | ||
lithomas1 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
minp : int, default 1 | ||
Minimum number of observations required per pair of columns | ||
to have a valid result. | ||
|
||
Returns | ||
------- | ||
numpy.ndarray | ||
lithomas1 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
Correlation matrix | ||
""" | ||
cdef: | ||
Py_ssize_t i, j, k, xi, yi, N, K | ||
ndarray[float64_t, ndim=2] result | ||
ndarray[float64_t, ndim=2] ranked_mat | ||
ndarray[uint8_t, ndim=2] mask | ||
float64_t currj | ||
ndarray[uint8_t, ndim=1] valid | ||
ndarray[int64_t] sorted_idxs | ||
ndarray[float64_t, ndim=1] col | ||
int64_t n_concordant | ||
int64_t total_concordant = 0 | ||
int64_t total_discordant = 0 | ||
float64_t kendall_tau | ||
int64_t n_obs | ||
const int64_t[:] labels_n | ||
|
||
N, K = (<object>mat).shape | ||
|
||
result = np.empty((K, K), dtype=np.float64) | ||
mask = np.isfinite(mat) | ||
|
||
ranked_mat = np.empty((N, K), dtype=np.float64) | ||
# For compatibility when calling rank_1d | ||
labels_n = np.zeros(N, dtype=np.int64) | ||
|
||
for i in range(K): | ||
ranked_mat[:, i] = rank_1d(mat[:, i], labels_n) | ||
|
||
for xi in range(K): | ||
sorted_idxs = ranked_mat[:, xi].argsort() | ||
ranked_mat = ranked_mat[sorted_idxs] | ||
mask = mask[sorted_idxs] | ||
for yi in range(xi + 1, K): | ||
valid = mask[:, xi] & mask[:, yi] | ||
if valid.sum() < minp: | ||
result[xi, yi] = NaN | ||
result[yi, xi] = NaN | ||
else: | ||
# Get columns and order second column using 1st column ranks | ||
if not valid.all(): | ||
col = ranked_mat[valid.nonzero()][:, yi] | ||
else: | ||
col = ranked_mat[:, yi] | ||
n_obs = col.shape[0] | ||
total_concordant = 0 | ||
total_discordant = 0 | ||
for j in range(n_obs - 1): | ||
currj = col[j] | ||
# Count num concordant and discordant pairs | ||
n_concordant = 0 | ||
for k in range(j, n_obs): | ||
if col[k] > currj: | ||
n_concordant += 1 | ||
total_concordant += n_concordant | ||
total_discordant += (n_obs-1-j-n_concordant) | ||
lithomas1 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
# Note: we do total_concordant+total_discordant here which is | ||
# equivalent to the C(n, 2), the total # of pairs, | ||
# listed on wikipedia | ||
kendall_tau = (total_concordant - total_discordant) / \ | ||
(total_concordant + total_discordant) | ||
result[xi, yi] = kendall_tau | ||
result[yi, xi] = kendall_tau | ||
|
||
if mask[:, xi].sum() > minp: | ||
result[xi, xi] = 1 | ||
else: | ||
result[xi, xi] = NaN | ||
|
||
return result | ||
|
||
|
||
# ---------------------------------------------------------------------- | ||
|
||
ctypedef fused algos_t: | ||
|
Uh oh!
There was an error while loading. Please reload this page.