1

嗨，我想计算相关性（r）值的 3D 矩阵的 p 值并计算完整 p 矩阵的 FDR 阈值（仅根据上三角 p 值计算 FDR）。3D 矩阵很大（208、3544、3544）。目前我使用这些 cython 函数来做到这一点。但是，我想让它更快，如果可能的话，让它使用更少的内存。你能帮我吗？

import numpy
cimport numpy
from scipy.special import stdtr 
DTYPE = numpy.float32
ctypedef numpy.float32_t DTYPE_t

def get_p_matrix_from_r_matrix_3D(DTYPE_t[:,:,::1] r_matrix_3D, int n):
    cdef unsigned int slice
    cdef numpy.ndarray[DTYPE_t, ndim=2] t_den
    cdef numpy.ndarray[DTYPE_t, ndim=2] t_
    cdef numpy.ndarray[DTYPE_t, ndim=2] r_matrix
    cdef numpy.ndarray[DTYPE_t, ndim=3] final_p = numpy.empty_like(r_matrix_3D)

    for slice in range(r_matrix_3D.shape[0]):
        r_matrix = numpy.asarray(r_matrix_3D[slice, :, :])
        t_den = numpy.sqrt((1 - r_matrix ** 2)/(n-2))
        t_ = r_matrix / t_den
        final_p[slice, :, :] = (1-stdtr(n - 2, numpy.abs(t_))) * 2

    return final_p

def get_FDG_threshold_man(numpy.ndarray[DTYPE_t, ndim=3] p_matrix, DTYPE_t alpha):
    index_array = numpy.ones((p_matrix.shape[0], p_matrix.shape[1], p_matrix.shape[2]), dtype=bool)
    cdef numpy.ndarray[DTYPE_t, ndim=1] p_s = p_matrix[numpy.triu(index_array, 1)]
    p_s.sort()
    cdef int m = len(p_s)
    cdef int k=0
    for i in reversed(range(m)):
        if (p_s[i] <= i/m*alpha) and (p_s[i+1] > (i+1)/m*alpha):
            k=i
            break
    if k != 0:
        p_th = p_s[k]
    else:
        p_th = 0
    print(p_th)
    return p_th
4

0 回答 0