0

我需要在 Julia 中做一些多项式回归。在 RI 中得到以下结果：

library(nnet)
data <- read.table("Dropbox/scripts/timeseries.txt",header=TRUE)
multinom(y~X1+X2,data)

# weights:  12 (6 variable)
initial  value 10985.024274 
iter  10 value 10438.503738
final  value 10438.503529 
converged
Call:
multinom(formula = y ~ X1 + X2, data = data)

Coefficients:
  (Intercept)        X1        X2
2   0.4877087 0.2588725 0.2762119
3   0.4421524 0.5305649 0.3895339

Residual Deviance: 20877.01 
AIC: 20889.01

这是我的数据

我的第一次尝试是使用Regression.jl。这个包的文档非常稀少，所以我不确定哪个类别用作基线，结果输出对应的参数等等。我在这里提出了一个问题来询问这些事情。

using DataFrames
using Regression
import Regression: solve, Options, predict
dat = readtable("timeseries.txt", separator='\t')
X = convert(Matrix{Float64},dat[:,2:3])
y = convert(Vector{Int64},dat[:,1])

ret = solve(mlogisticreg(X',y,3), reg=ZeroReg(),  options=Options(verbosity=:iter))

结果是

julia> ret.sol
3x2 Array{Float64,2}:
 -0.573027  -0.531819
  0.173453   0.232029
  0.399575   0.29979

但同样，我不确定这对应于什么。

接下来我尝试了 Python 的 SciKitLearn 的 Julia 包装器：

using ScikitLearn
@sk_import linear_model: LogisticRegression

model = ScikitLearn.fit!(LogisticRegression(multi_class="multinomial", solver = "lbfgs"), X, y)
model[:coef_]
3x2 Array{Float64,2}:
 -0.261902    -0.220771 
 -0.00453731   0.0540354
  0.266439     0.166735

但我还没有弄清楚如何从这个模型中提取系数。用系数更新。这些看起来也不像 R 结果。

任何试图复制 R 的结果的帮助将不胜感激（使用任何包！）

请注意，响应变量只是离散化的时滞响应，即

julia> dat[1:3,:]
3x3 DataFrames.DataFrame
| Row | y | X1 | X2 |
|-----|---|----|----|
| 1   | 3 | 1  | 0  |
| 2   | 3 | 0  | 1  |
| 3   | 1 | 0  | 1  |

对于第 2 行，您可以看到响应 (0, 1) 表示之前的观察结果是 3。类似地，(1,0) 表示之前的观察结果是 2，(0,0) 表示之前的观察结果是 1。

更新：对于 Regression.jl，默认情况下它似乎不适合截距（他们称其为“偏差”而不是截距）。通过添加这个术语，我们得到与 python 非常相似的结果（虽然不确定第三列是什么......）

julia> ret = solve(mlogisticreg(X',y,3, bias=1.0), reg=ZeroReg(),  options=Options(verbosity=:iter))
julia> ret.sol
3x3 Array{Float64,2}:
 -0.263149    -0.221923   -0.309949
 -0.00427033   0.0543008   0.177753
  0.267419     0.167622    0.132196

更新：由于模型系数无法识别，我不应该期望它们在这些不同的实现中是相同的。但是，预测的概率应该是相同的，实际上它们是相同的（使用 R、Regression.jl 或 ScikitLearn）。

4

0 回答 0