小编Jos*_*sef的帖子

我是Python的新手,来自R世界.我正在尝试使用SciPy将分布拟合到样本数据并取得了很好的成功.我可以distribution.fit(data)回报理智的结果.我一直无法做的是创建拟合优度统计数据,这是我习惯使用fitdistrplusR中的包.是否有一种通用的方法来比较SciPy中许多不同发行版的"最佳拟合"？

我正在寻找像Kolmogorov-Smirnov测试或Cramer-von Mises或Anderson-darling测试的东西

我想学习使用Python的statsmodels库中的普通最小二乘模型,描述在这里.

sm.OLS.fit()返回学习的模型.有没有办法将其保存到文件并重新加载？我的训练数据很大,学习模型大约需要半分钟.所以我想知道OLS模型中是否存在任何保存/加载功能.

我repr()在模型对象上尝试了该方法,但它没有返回任何有用的信息.

我正在尝试对德国信用数据运行logit回归(www4.stat.ncsu.edu/~boos/var.select/german.credit.html).为了测试代码,我只使用了数值变量,并尝试使用以下代码对结果进行回归.

import pandas as pd
import statsmodels.api as sm
import pylab as pl
import numpy as np

df = pd.read_csv("germandata.txt",delimiter=' ')
df.columns = ["chk_acc","duration","history","purpose","amount","savings_acc","employ_since","install_rate","pers_status","debtors","residence_since","property","age","other_plans","housing","existing_credit","job","no_people_liab","telephone","foreign_worker","admit"]

#pls note that I am only retaining numeric variables
cols_to_keep = ['admit','duration', 'amount', 'install_rate','residence_since','age','existing_credit','no_people_liab']

# rank of cols_to_keep is 8
print np.linalg.matrix_rank(df[cols_to_keep].values)
data = df[cols_to_keep]

data['intercept'] = 1.0

train_cols = data.columns[1:]

#to check the rank of train_cols, which in this case is 8
print np.linalg.matrix_rank(data[train_cols].values)

#fit logit model
logit = sm.Logit(data['admit'], data[train_cols])
result = logit.fit()

检查数据时,所有8.0列都是独立的.尽管如此,我得到奇异矩阵错误.你能帮忙吗？ …

我有一个模型,定义如下:

import statsmodels.formula.api as smf
model = smf.glm(formula="A ~ B + C + D", data=data, family=sm.families.Poisson()).fit()

该模型具有如下所示的系数:

Intercept   0.319813
C[T.foo]   -1.058058
C[T.bar]   -0.749859
D[T.foo]    0.217136
D[T.bar]    0.404791
B           0.262614

我可以抓住的值Intercept,并B通过做model.params.Intercept和model.params.B,但我不能让每个值C和D.

我试过model.params.C[T.foo],例如,我得到了错误.

我如何从模型中获取特定值？

我试图为Genz和Bretz之后的多变量t分布编程cdf的算法,R中的参考包是mvtnorm.

当我测试我的功能时,我发现我的数字不匹配.在以下示例中,从mvtnorm帮助调整,多变量t随机变量具有独立的组件.所以积分应该只是3个独立概率的乘积

> lower <- -1
> upper <- 3
> df <- 4
> corr <- diag(3)
> delta <- rep(0, 3)
> pmvt(lower=lower, upper=upper, delta=delta, df=df, corr=corr)
[1] 0.5300413
attr(,"error")
[1] 4.321136e-05
attr(,"msg")
[1] "Normal Completion"

报告的误差是4e-5,该误差与独立概率的乘积相比较

> (pt(upper, df) - pt(lower, df))**3
[1] 0.4988254

是

0.5300413 - 0.4988254 = 0.0312159

与R mvtnorm相比,我在自己的代码中得到的差异在于大约相同范围内的各种示例.

我大部分都是R的初学者.那么,我做错了什么或出了什么问题？

(我没有在R-help邮件列表上注册,所以我在这里试试.)

更新:正如pchalasani解释的那样,我的统计数据是错误的,我自己的代码中的错误是在一些辅助函数中,而不是在t分发代码中.看到不相关的好方法并不意味着独立,正在考虑条件分布.以下是四分位数的独立双变量随机变量(10000个样本)的列频率%*100(以列变量为条件的分布).

双变量不相关的正常变量

([[26, 25, 24, 23],
  [24, 23, 24, 25],
  [24, 27, 24, 24],
  [24, 23, 26, 25]])

双变量不相关的变量

([[29, 20, 22, 29], …

使用 statesmodels 的逻辑回归模型：

log_reg = st.logit(formula = 'label ~ pregnant + glucose + bp + insulin + bmi + pedigree + age', data=pima).fit()

有没有简单的方法来写公式的第二部分（怀孕+葡萄糖+血压+胰岛素+体重指数+血统+年龄）？这里必须明确提及所有列。如果超过100列，就很难写，而且语句也会很长。

我想要一个与之相关的系数和Newey-West标准误.

我正在寻找Python库(理想情况下,但任何工作解决方案都很好)可以做以下R代码正在做的事情:

library(sandwich)
library(lmtest)

a <- matrix(c(1,3,5,7,4,5,6,4,7,8,9))
b <- matrix(c(3,5,6,2,4,6,7,8,7,8,9))

temp.lm = lm(a ~ b)

temp.summ <- summary(temp.lm)
temp.summ$coefficients <- unclass(coeftest(temp.lm, vcov. = NeweyWest))

print (temp.summ$coefficients)

结果:

             Estimate Std. Error   t value  Pr(>|t|)
(Intercept) 2.0576208  2.5230532 0.8155281 0.4358205
b           0.5594796  0.4071834 1.3740235 0.2026817

我得到系数并与它们相关的标准误差.

我看到statsmodels.stats.sandwich_covariance.cov_hac模块,但我不知道如何使它与OLS一起工作.

我一直在尝试复制这个矩阵（R )Julia

data = read.csv("https://gist.githubusercontent.com/TJhon/158daa0c2dd06010d01a72dae2af8314/raw/61df065c98ec90b9ea3b8598d1996fb5371a64aa/rnd.csv")\n\nhead(model.matrix(y ~ (x1 + x2 + x3 + x4 + x5)^2, data), 3)\n#>   (Intercept)         x1         x2         x3        x4 x5       x1:x2\n#> 1           1 -0.3007225 -1.3710894  0.3423409  1.322547  2  0.41231744\n#> 2           1  0.4674170  0.8728939  0.9534157 -1.007083  1  0.40800548\n#> 3           1  0.2085316 -0.3657995 -0.3043694 -1.036938  4 -0.07628076\n#>         x1:x3      x1:x4      x1:x5      x2:x3      x2:x4      x2:x5      x3:x4\n#> 1 -0.10294961 -0.3977198 -0.6014450 -0.4693799 -1.8133307 -2.7421787  0.4527620\n#> 2  0.44564276 -0.4707279  0.4674170  0.8322308 -0.8790769  0.8728939 -0.9601690\n#> 3 -0.06347064 -0.2162343  0.8341265  0.1113382  0.3793113 …

所以我尝试使用 python 的 statsmodels.api 对二进制结果进行逻辑回归进行预测。我按照教程使用 Logit。当我尝试对测试数据集进行预测时，每个记录的输出都是 0 到 1 之间的小数。它不应该给我零和一吗？或者我是否必须使用圆形函数或其他东西来转换它们。

请原谅这个问题的新手。我正在凝视我的旅程。

我有一些数据集：titanic

在R中这样做

glm(Survived ~ Sex, titanic, family = "binomial")

我懂了

(Intercept)     SexMale 
1.124321   -2.477825

R以幸存为阳性。

但是当我在Python中做同样的事情时

sm.formula.glm("Survived ~ Sex", family=sm.families.Binomial(), data=titanic).fit()

我得到负面结果：即Python 不能幸免为正面结果。

如何调整Python的glm函数行为，使其返回与R相同的结果？