def gp_predict(gpmodel, para_array_rescaled):
    """Returns a Gaussian Process prediction of PCA weights

    Parameters
    ----------

    gpmodel: GPflow object
        Trained GPflow model corresponding to a particular redshift

    para_array_rescaled: ndarray of shape (5, )
        Rescaled Cosmological parameters


    Returns
    _______

    W_predArray: float or ndarray of shape (nRankMax, )
        Mean GP prediction of the PCA weights

    W_varArray: float or ndarray of shape (nRankMax, )
        Variance in the GP prediction of the PCA weights

    """

    m1p = gpmodel.predict_f(para_array_rescaled)

    # [0] is the mean and [1] is the std of the prediction
    W_predArray = m1p[0]
    W_varArray = m1p[1]

    return W_predArray, W_varArray