在多变量回归的情况下，@JARH的答案可能有错误。 (我没有足够的声誉来评论。)

在下面一行:

p_values = [2 * (1-stats.t.cdf (np.abs(我),(len (newX) 1)))我在ts_b),

t值遵循degree len(newX)-1的卡方分布，而不是遵循degree len(newX)-len(newX.columns)-1的卡方分布。

所以这应该是:

p_values = [2 * (1-stats.t.cdf (np.abs(我),(len (newX) len (newX.columns) 1)))我在ts_b)

(详见OLS回归的t值)

获取p值的一个简单方法是使用statmodels回归:

import statsmodels.api as sm
mod = sm.OLS(Y,X)
fii = mod.fit()
p_values = fii.summary2().tables[1]['P>|t|']

你可以得到一系列你可以操作的p值(例如，通过计算每个p值来选择你想要保持的顺序):

scikit-learn的线性回归不计算这个信息，但你可以很容易地扩展类来做:

from sklearn import linear_model
from scipy import stats
import numpy as np


class LinearRegression(linear_model.LinearRegression):
    """
    LinearRegression class after sklearn's, but calculate t-statistics
    and p-values for model coefficients (betas).
    Additional attributes available after .fit()
    are `t` and `p` which are of the shape (y.shape[1], X.shape[1])
    which is (n_features, n_coefs)
    This class sets the intercept to 0 by default, since usually we include it
    in X.
    """

    def __init__(self, *args, **kwargs):
        if not "fit_intercept" in kwargs:
            kwargs['fit_intercept'] = False
        super(LinearRegression, self)\
                .__init__(*args, **kwargs)

    def fit(self, X, y, n_jobs=1):
        self = super(LinearRegression, self).fit(X, y, n_jobs)

        sse = np.sum((self.predict(X) - y) ** 2, axis=0) / float(X.shape[0] - X.shape[1])
        se = np.array([
            np.sqrt(np.diagonal(sse[i] * np.linalg.inv(np.dot(X.T, X))))
                                                    for i in range(sse.shape[0])
                    ])

        self.t = self.coef_ / se
        self.p = 2 * (1 - stats.t.cdf(np.abs(self.t), y.shape[0] - X.shape[1]))
        return self

从这里偷来的。

您应该看一看Python中用于这种统计分析的统计模型。

稍微了解一下线性回归理论，下面是我们计算系数估计器(随机变量)的p值所需的总结，以检查它们是否显著(通过拒绝相应的零假设):

现在，让我们用下面的代码段计算p值:

import numpy as np 
# generate some data 
np.random.seed(1)
n = 100
X = np.random.random((n,2))
beta = np.array([-1, 2])
noise = np.random.normal(loc=0, scale=2, size=n)
y = X@beta + noise

用scikit-learn从上面的公式计算p值:

# use scikit-learn's linear regression model to obtain the coefficient estimates
from sklearn.linear_model import LinearRegression
reg = LinearRegression().fit(X, y)
beta_hat = [reg.intercept_] + reg.coef_.tolist()
beta_hat
# [0.18444290873001834, -1.5879784718284842, 2.5252138207251904]

# compute the p-values
from scipy.stats import t
# add ones column
X1 = np.column_stack((np.ones(n), X))
# standard deviation of the noise.
sigma_hat = np.sqrt(np.sum(np.square(y - X1@beta_hat)) / (n - X1.shape[1]))
# estimate the covariance matrix for beta 
beta_cov = np.linalg.inv(X1.T@X1)
# the t-test statistic for each variable from the formula from above figure
t_vals = beta_hat / (sigma_hat * np.sqrt(np.diagonal(beta_cov)))
# compute 2-sided p-values.
p_vals = t.sf(np.abs(t_vals), n-X1.shape[1])*2 
t_vals
# array([ 0.37424023, -2.36373529,  3.57930174])
p_vals
# array([7.09042437e-01, 2.00854025e-02, 5.40073114e-04])

使用statmodels计算p值:

import statsmodels.api as sm
X1 = sm.add_constant(X)
model = sm.OLS(y, X2)
model = model.fit()
model.tvalues
# array([ 0.37424023, -2.36373529,  3.57930174])
# compute p-values
t.sf(np.abs(model.tvalues), n-X1.shape[1])*2 
# array([7.09042437e-01, 2.00854025e-02, 5.40073114e-04])  

model.summary()

从上面可以看出，两种情况下计算的p值完全相同。

在多变量回归的情况下，@JARH的答案可能有错误。 (我没有足够的声誉来评论。)

在下面一行:

p_values = [2 * (1-stats.t.cdf (np.abs(我),(len (newX) 1)))我在ts_b),

t值遵循degree len(newX)-1的卡方分布，而不是遵循degree len(newX)-len(newX.columns)-1的卡方分布。

所以这应该是:

p_values = [2 * (1-stats.t.cdf (np.abs(我),(len (newX) len (newX.columns) 1)))我在ts_b)

(详见OLS回归的t值)

另外一个已经提出的选择是使用排列测试。将y的值洗牌后对模型进行N次拟合，计算拟合后模型的系数相对于原模型的系数值较大(单侧检验)或绝对值较大(双面检验)的比例。这些比例就是p值。