wedin_bound_svec_resampling.py 文件源码

def norms_svec_resampling(X, orthogonal_basis, rank, num_samples=1000):
    """
    Resampling procedure described in AJIVE paper for Wedin bound

    Parameters
    ---------
    orthogonal_basis: basis vectors for the orthogonal complement
    of the score space

    X: the observed data

    rank: number of columns to resample

    num_samples: how many resamples

    Output
    ------
    an array of the resampled norms
    """
    # TODO: rename some arguments e.g. orthogonal_basis -- this is not really
    # a basis, also resampled_projection

    resampled_norms = [0]*num_samples

    for s in range(num_samples):

        # sample columns from orthogonal basis
        sampled_col_index = np.random.choice(orthogonal_basis.shape[1], size=rank, replace=True)
        resampled_basis = orthogonal_basis[:, sampled_col_index]
        # ^ this is V* from AJIVE p12

        # project observed data

        # resampled_projection = np.dot(X, resampled_basis)
        resampled_projection = safe_sparse_dot(X, resampled_basis)


        # compute resampled operator L2 norm
        resampled_norms[s] = np.linalg.norm(resampled_projection, ord=2)

    return resampled_norms