• Archive by category "R语言实践"

Blog Archives

R语言解读多元线性回归模型

R的极客理想系列文章,涵盖了R的思想,使用,工具,创新等的一系列要点,以我个人的学习和体验去诠释R的强大。

R语言作为统计学一门语言,一直在小众领域闪耀着光芒。直到大数据的爆发,R语言变成了一门炙手可热的数据分析的利器。随着越来越多的工程背景的人的加入,R语言的社区在迅速扩大成长。现在已不仅仅是统计领域,教育,银行,电商,互联网….都在使用R语言。

要成为有理想的极客,我们不能停留在语法上,要掌握牢固的数学,概率,统计知识,同时还要有创新精神,把R语言发挥到各个领域。让我们一起动起来吧,开始R的极客理想。

关于作者:

  • 张丹(Conan), 程序员R,Nodejs,Java
  • weibo:@Conan_Z
  • blog: http://blog.fens.me
  • email: bsspirit@gmail.com

转载请注明出处:
http://blog.fens.me/r-multi-linear-regression

reg-multi-liner

前言

本文接上一篇R语言解读一元线性回归模型。在许多生活和工作的实际问题中,影响因变量的因素可能不止一个,比如对于知识水平越高的人,收入水平也越高,这样的一个结论。这其中可能包括了因为更好的家庭条件,所以有了更好的教育;因为在一线城市发展,所以有了更好的工作机会;所处的行业赶上了大的经济上行周期等。要想解读这些规律,是复杂的、多维度的,多元回归分析方法更适合解读生活的规律。

由于本文为非统计的专业文章,所以当出现与教课书不符的描述,请以教课书为准。本文力求用简化的语言,来介绍多元线性回归的知识,同时配合R语言的实现。

目录

  1. 多元线性回归介绍
  2. 元线性回归建模
  3. 模型优化
  4. 案例:黑色系期货日K线数据验证

1. 多元线性回归介绍

对比一元线性回归,多元线性回归是用来确定2个或2个以上变量间关系的统计分析方法。多元线性回归的基本的分析方法与一元线性回归方法是类似的,我们首先需要对选取多元数据集并定义数学模型,然后进行参数估计,对估计出来的参数进行显著性检验,残差分析,异常点检测,最后确定回归方程进行模型预测。

由于多元回归方程有多个自变量,区别于一元回归方程,有一项很重要的操作就是自变量的优化,挑选出相关性最显著的自变量,同时去除不显著的自变量。在R语言中,有很方便地用于优化函数,可以很好的帮助我们来改进回归模型。

下面就开始多元线性回归的建模过程。

2. 多元线性回归建模

做过商品期货研究的人,都知道黑色系品种是具有产业链上下游的关系。铁矿石是炼钢的原材料,焦煤和焦炭是炼钢的能源资源,热卷即热轧卷板是以板坯为原料经加热后制成的钢板,螺纹钢是表面带肋的钢筋。

由于有产业链的关系,假设我们想要预测螺纹钢的价格,那么影响螺纹钢价格的因素可以会涉及到原材料,能源资源和同类材料等。比如,铁矿石价格如果上涨,螺纹钢就应该要涨价了。

2.1 数据集和数学模型

先从数据开始介绍,这次的数据集,我选择的期货黑色系的品种的商品期货,包括了大连期货交易所的 焦煤(JM),焦炭(J),铁矿石(I),上海期货交易所的 螺纹钢(RU) 和 热卷(HC)。

数据集为2016年3月15日,当日白天开盘的交易数据,为黑色系的5个期货合约的分钟线的价格数据。


# 数据集已存在df变量中
> head(df,20)
                       x1    x2    x3   x4    y
2016-03-15 09:01:00 754.5 616.5 426.5 2215 2055
2016-03-15 09:02:00 752.5 614.5 423.5 2206 2048
2016-03-15 09:03:00 753.0 614.0 423.0 2199 2044
2016-03-15 09:04:00 752.5 613.0 422.5 2197 2040
2016-03-15 09:05:00 753.0 615.5 424.0 2198 2043
2016-03-15 09:06:00 752.5 614.5 422.0 2195 2040
2016-03-15 09:07:00 752.0 614.0 421.5 2193 2036
2016-03-15 09:08:00 753.0 615.0 422.5 2197 2043
2016-03-15 09:09:00 754.0 615.5 422.5 2197 2041
2016-03-15 09:10:00 754.5 615.5 423.0 2200 2044
2016-03-15 09:11:00 757.0 616.5 423.0 2201 2045
2016-03-15 09:12:00 756.0 615.5 423.0 2200 2044
2016-03-15 09:13:00 755.5 615.0 423.0 2197 2042
2016-03-15 09:14:00 755.5 615.0 423.0 2196 2042
2016-03-15 09:15:00 756.0 616.0 423.5 2200 2045
2016-03-15 09:16:00 757.5 616.0 424.0 2205 2052
2016-03-15 09:17:00 758.5 618.0 424.0 2204 2051
2016-03-15 09:18:00 759.5 618.5 424.0 2205 2053
2016-03-15 09:19:00 759.5 617.5 424.5 2206 2053
2016-03-15 09:20:00 758.5 617.5 423.5 2201 2050

数据集包括有6列:

  • 索引, 为时间
  • x1, 为焦炭(j1605)合约的1分钟线的报价数据
  • x2, 为焦煤(jm1605)合约的1分钟线的报价数据
  • x3, 为铁矿石(i1605)合约的1分钟线的报价数据
  • x4, 为热卷(hc1605)合约的1分钟线的报价数据
  • y, 为螺纹钢(rb1605)合约的1分钟线的报价数据

假设螺纹钢的价格与其他4个商品的价格有线性关系,那么我们建立以螺纹钢为因变量,以焦煤、焦炭、铁矿石和热卷的为自变量的多元线性回归模型。用公式表示为:

y = a + b * x1 + c * x2 + d * x3 + e * x4 + ε
  • y,为因变量,螺纹钢
  • x1,为自变量,焦煤
  • x2,为自变量,焦炭
  • x3,为自变量,铁矿石
  • x4,为自变量,热卷
  • a,为截距
  • b,c,d,e,为自变量系数
  • ε, 为残差,是其他一切不确定因素影响的总和,其值不可观测。假定ε服从正态分布N(0,σ^2)。

通过对多元线性回归模型的数学定义,接下来让我们利用数据集做多元回归模型的参数估计。

2.2. 回归参数估计

上面公式中,回归参数 a, b, c, d,e都是我们不知道的,参数估计就是通过数据来估计出这些参数,从而确定自变量和因变量之前的关系。我们的目标是要计算出一条直线,使直线上每个点的Y值和实际数据的Y值之差的平方和最小,即(Y1实际-Y1预测)^2+(Y2实际-Y2预测)^2+ …… +(Yn实际-Yn预测)^2 的值最小。参数估计时,我们只考虑Y随X自变量的线性变化的部分,而残差ε是不可观测的,参数估计法并不需要考虑残差。

类似于一元线性回归,我们用R语言来实现对数据的回归模型的参数估计,用lm()函数来实现多元线性回归的建模过程。


# 建立多元线性回归模型
> lm1<-lm(y~x1+x2+x3+x4,data=df)

# 打印参数估计的结果
> lm1

Call:
lm(formula = y ~ x1 + x2 + x3 + x4, data = df)

Coefficients:
(Intercept)           x1           x2           x3           x4  
   212.8780       0.8542       0.6672      -0.6674       0.4821  

这样我们就得到了y和x关系的方程。

y = 212.8780 + 0.8542 * x1 + 0.6672 * x2 - 0.6674 * x3 + 0.4821 * x4

2.3. 回归方程的显著性检验

参考一元线性回归的显著性检验,多元线性回归的显著性检验,同样是需要经过 T检验,F检验,和R^2(R平方)相关系统检验。在R语言中这三种检验的方法都已被实现,我们只需要把结果解读,我们可以summary()函数来提取模型的计算结果。


> summary(lm1)

Call:
lm(formula = y ~ x1 + x2 + x3 + x4, data = df)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.9648 -1.3241 -0.0319  1.2403  5.4194 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 212.87796   58.26788   3.653 0.000323 ***
x1            0.85423    0.10958   7.795 2.50e-13 ***
x2            0.66724    0.12938   5.157 5.57e-07 ***
x3           -0.66741    0.15421  -4.328 2.28e-05 ***
x4            0.48214    0.01959  24.609  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.028 on 221 degrees of freedom
Multiple R-squared:  0.9725,	Adjusted R-squared:  0.972 
F-statistic:  1956 on 4 and 221 DF,  p-value: < 2.2e-16
  • T检验:所自变量都是非常显著***
  • F检验:同样是非常显著,p-value < 2.2e-16
  • 调整后的R^2:相关性非常强为0.972

最后,我们通过的回归参数的检验与回归方程的检验,得到最后多元线性回归方程为:


y = 212.87796 + 0.85423 * x1 + 0.66724 * x2 - 0.66741 * x3 + 0.48214 * x4

即

螺纹钢 = 212.87796 + 0.85423 * 焦炭 + 0.66724 * 焦煤 - 0.66741 * 铁矿石 + 0.48214 * 热卷

2.4 残差分析和异常点检测

在得到的回归模型进行显著性检验后,还要在做残差分析(预测值和实际值之间的差),检验模型的正确性,残差必须服从正态分布N(0,σ^2)。直接用plot()函数生成4种用于模型诊断的图形,进行直观地分析。


> par(mfrow=c(2,2))
> plot(lm1)

m01

  • 残差和拟合值(左上),残差和拟合值之间数据点均匀分布在y=0两侧,呈现出随机的分布,红色线呈现出一条平稳的曲线并没有明显的形状特征。
  • 残差QQ图(右上),数据点按对角直线排列,趋于一条直线,并被对角直接穿过,直观上符合正态分布。
  • 标准化残差平方根和拟合值(左下),数据点均匀分布在y=0两侧,呈现出随机的分布,红色线呈现出一条平稳的曲线并没有明显的形状特征。
  • 标准化残差和杠杆值(右下),没有出现红色的等高线,则说明数据中没有特别影响回归结果的异常点。

结论,没有明显的异常点,残差符合假设条件。

2.5. 模型预测

我们得到了多元线性回归方程的公式,就可以对数据进行预测了。我们可以用R语言的predict()函数来计算预测值y0和相应的预测区间,并把实际值和预测值一起可视化化展示。


> par(mfrow=c(1,1))  #设置画面布局

# 预测计算
> dfp<-predict(lm1,interval="prediction")

# 打印预测时
> head(dfp,10)
                fit      lwr      upr
2014-03-21 3160.526 3046.425 3274.626
2014-03-24 3193.253 3078.868 3307.637
2014-03-25 3240.389 3126.171 3354.607
2014-03-26 3228.565 3114.420 3342.710
2014-03-27 3222.528 3108.342 3336.713
2014-03-28 3262.399 3148.132 3376.666
2014-03-31 3291.996 3177.648 3406.344
2014-04-01 3305.870 3191.447 3420.294
2014-04-02 3275.370 3161.018 3389.723
2014-04-03 3297.358 3182.960 3411.755

# 合并数据
> mdf<-merge(df$y,dfp)	 

# 画图
> draw(mdf)

m02

图例说明

  • y, 实际价格,红色线
  • fit, 预测价格,绿色线
  • lwr,预测最低价,蓝色线
  • upr,预测最高价,紫色线

从图中看出,实际价格y和预测价格fit,在大多数的时候都是很贴近的。我们的一个模型就训练好了!

3. 模型优化

上文中,我们已经很顺利的发现了一个非常不错的模型。如果要进行模型优化,可以用R语言中update()函数进行模型的调整。我们首先检查一下每个自变量x1,x2,x3,x4和因变量y之间的关系。

pairs(as.data.frame(df))

m03

从图中,我们可以发现x2与Y的关系,可能是最偏离线性的。那么,我们尝试对多元线性回归模型进行调整,从原模型中去掉x2变量。



# 模型调整
> lm2<-update(lm1, .~. -x2)

> summary(lm2)

Call:
lm(formula = y ~ x1 + x3 + x4, data = df)

Residuals:
    Min      1Q  Median      3Q     Max 
-6.0039 -1.3842  0.0177  1.3513  4.8028 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 462.47104   34.26636   13.50  < 2e-16 ***
x1            1.08728    0.10543   10.31  < 2e-16 ***
x3           -0.40788    0.15394   -2.65  0.00864 ** 
x4            0.42582    0.01718   24.79  < 2e-16 ***
---
Signif. codes:  
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.142 on 222 degrees of freedom
Multiple R-squared:  0.9692,	Adjusted R-squared:  0.9688 
F-statistic:  2330 on 3 and 222 DF,  p-value: < 2.2e-16

当把自变量x2去掉后,自变量x3的T检验反而变大了,同时Adjusted R-squared变小了,所以我们这次调整是有问题的。

如果通过生产和原材料的内在逻辑分析,焦煤与焦炭属于上下游关系。焦煤是生产焦炭的一种原材料,焦炭是焦煤与其他炼焦煤经过配煤焦化形成的产品,一般生产 1 吨焦炭需要1.33 吨炼焦煤,其中焦煤至少占 30% 。

我们把焦煤 和 焦炭的关系改变一下,增加x1*x2的关系匹配到模型,看看效果。


# 模型调整
> lm3<-update(lm1, .~. + x1*x2)
> summary(lm3)

Call:
lm(formula = y ~ x1 + x2 + x3 + x4 + x1:x2, data = df)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.8110 -1.3501 -0.0595  1.2019  5.3884 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) 7160.32231 7814.50048   0.916    0.361    
x1            -8.45530   10.47167  -0.807    0.420    
x2           -10.58406   12.65579  -0.836    0.404    
x3            -0.64344    0.15662  -4.108 5.63e-05 ***
x4             0.48363    0.01967  24.584  < 2e-16 ***
x1:x2          0.01505    0.01693   0.889    0.375    
---
Signif. codes:  
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.029 on 220 degrees of freedom
Multiple R-squared:  0.9726,	Adjusted R-squared:  0.972 
F-statistic:  1563 on 5 and 220 DF,  p-value: < 2.2e-16

从结果中发现,增加了x1*x2列后,原来的x1,x2和Intercept的T检验都不显著。继续调整模型,从模型中去掉x1,x2两个自变量。


# 模型调整
> lm4<-update(lm3, .~. -x1-x2)
> summary(lm4)

Call:
lm(formula = y ~ x3 + x4 + x1:x2, data = df)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.9027 -1.2516 -0.0167  1.2748  5.8683 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  6.950e+02  1.609e+01  43.183  < 2e-16 ***
x3          -6.284e-01  1.530e-01  -4.108 5.61e-05 ***
x4           4.959e-01  1.785e-02  27.783  < 2e-16 ***
x1:x2        1.133e-03  9.524e-05  11.897  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.035 on 222 degrees of freedom
Multiple R-squared:  0.9722,	Adjusted R-squared:  0.9718 
F-statistic:  2588 on 3 and 222 DF,  p-value: < 2.2e-16

从调整后的结果来看,效果还不错。不过,也并没有比最初的模型有所提高。

对于模型调整的过程,如果我们手动调整测试时,一般都会基于业务知识来操作。如果是按照数据指标来计算,我们可以用R语言中提供的逐步回归的优化方法,通过AIC指标来判断是否需要参数优化。


#对lm1模型做逐步回归
> step(lm1)
Start:  AIC=324.51
y ~ x1 + x2 + x3 + x4

       Df Sum of Sq    RSS    AIC
               908.8 324.51
- x3    1     77.03  985.9 340.90
- x2    1    109.37 1018.2 348.19
- x1    1    249.90 1158.8 377.41
- x4    1   2490.56 3399.4 620.65

Call:
lm(formula = y ~ x1 + x2 + x3 + x4, data = df)

Coefficients:
(Intercept)           x1           x2           x3           x4  
   212.8780       0.8542       0.6672      -0.6674       0.4821  

通过计算AIC指标,lm1的模型AIC最小时为324.51,每次去掉一个自变量都会让AIC的值变大,所以我们还是不调整比较好。

对刚才的lm3模型做逐步回归的模型调整。


#对lm3模型做逐步回归
> step(lm3)
Start:  AIC=325.7               #当前AIC
y ~ x1 + x2 + x3 + x4 + x1:x2

        Df Sum of Sq    RSS    AIC
- x1:x2  1      3.25  908.8 324.51
                905.6 325.70
- x3     1     69.47  975.1 340.41
- x4     1   2487.86 3393.5 622.25

Step:  AIC=324.51               #去掉x1*x2项的AIC
y ~ x1 + x2 + x3 + x4

       Df Sum of Sq    RSS    AIC
               908.8 324.51
- x3    1     77.03  985.9 340.90
- x2    1    109.37 1018.2 348.19
- x1    1    249.90 1158.8 377.41
- x4    1   2490.56 3399.4 620.65

Call:
lm(formula = y ~ x1 + x2 + x3 + x4, data = df)

Coefficients:
(Intercept)           x1           x2           x3           x4  
   212.8780       0.8542       0.6672      -0.6674       0.4821  

通过AIC的判断,去掉X1*X2项后AIC最小,最后的检验结果告诉我们,还是原初的模型是最好的。

4. 案例:黑色系期货日K线数据验证

最后,我们用上面5个期货合约的日K线数据测试一下,找到多元回归关系。


> lm9<-lm(y~x1+x2+x3+x4,data=df)  # 日K线数据
> summary(lm9)

Call:
lm(formula = y ~ x1 + x2 + x3 + x4, data = df)

Residuals:
     Min       1Q   Median       3Q      Max 
-173.338  -37.470    3.465   32.158  178.982 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 386.33482   31.07729  12.431  < 2e-16 ***
x1            0.75871    0.07554  10.045  < 2e-16 ***
x2           -0.62907    0.14715  -4.275 2.24e-05 ***
x3            1.16070    0.05224  22.219  < 2e-16 ***
x4            0.46461    0.02168  21.427  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 57.78 on 565 degrees of freedom
Multiple R-squared:  0.9844,	Adjusted R-squared:  0.9843 
F-statistic:  8906 on 4 and 565 DF,  p-value: < 2.2e-16

数据集的基本统计信息。


> summary(df)
     Index                           x1               x2       
 Min.   :2014-03-21 00:00:00   Min.   : 606.5   Min.   :494.0  
 1st Qu.:2014-10-21 06:00:00   1st Qu.: 803.5   1st Qu.:613.1  
 Median :2015-05-20 12:00:00   Median : 939.0   Median :705.8  
 Mean   :2015-05-21 08:02:31   Mean   : 936.1   Mean   :695.3  
 3rd Qu.:2015-12-16 18:00:00   3rd Qu.:1075.0   3rd Qu.:773.0  
 Max.   :2016-07-25 00:00:00   Max.   :1280.0   Max.   :898.0  

       x3              x4             y       
 Min.   :284.0   Min.   :1691   Min.   :1626  
 1st Qu.:374.1   1st Qu.:2084   1st Qu.:2012  
 Median :434.0   Median :2503   Median :2378  
 Mean   :476.5   Mean   :2545   Mean   :2395  
 3rd Qu.:545.8   3rd Qu.:2916   3rd Qu.:2592  
 Max.   :825.0   Max.   :3480   Max.   :3414  

m04

对于日K线数据,黑色系的5个品种,同样具有非常强的相关关系,那么我们就可以把这个结论应用到实际的交易中了。

本文通过多元回归的统计分析方法,介绍多元回归在金融市场的基本应用。我们通过建立因变量和多个自变量的模型,从而发现生活中更复杂的规律,并建立有效的验证指标。让我们我们的技术优势,去金融市场抢钱吧。

转载请注明出处:
http://blog.fens.me/r-multi-linear-regression

打赏作者

R语言解读一元线性回归模型

R的极客理想系列文章,涵盖了R的思想,使用,工具,创新等的一系列要点,以我个人的学习和体验去诠释R的强大。

R语言作为统计学一门语言,一直在小众领域闪耀着光芒。直到大数据的爆发,R语言变成了一门炙手可热的数据分析的利器。随着越来越多的工程背景的人的加入,R语言的社区在迅速扩大成长。现在已不仅仅是统计领域,教育,银行,电商,互联网….都在使用R语言。

要成为有理想的极客,我们不能停留在语法上,要掌握牢固的数学,概率,统计知识,同时还要有创新精神,把R语言发挥到各个领域。让我们一起动起来吧,开始R的极客理想。

关于作者:

  • 张丹(Conan), 程序员R,Nodejs,Java
  • weibo:@Conan_Z
  • blog: http://blog.fens.me
  • email: bsspirit@gmail.com

转载请注明出处:
http://blog.fens.me/r-linear-regression/

reg-liner

前言

在我们的日常生活中,存在大量的具有相关性的事件,比如大气压和海拔高度,海拔越高大气压强越小;人的身高和体重,普遍来看越高的人体重也越重。还有一些可能存在相关性的事件,比如知识水平越高的人,收入水平越高;市场化的国家经济越好,则货币越强势,反而全球经济危机,黄金等避险资产越走强。

如果我们要研究这些事件,找到不同变量之间的关系,我们就会用到回归分析。一元线性回归分析是处理两个变量之间关系的最简单模型,是两个变量之间的线性相关关系。让我们一起发现生活中的规律吧。

由于本文为非统计的专业文章,所以当出现与教课书不符的描述,请以教课书为准。本文力求用简化的语言,来介绍一元线性回归的知识,同时配合R语言的实现。

目录

  1. 一元线性回归介绍
  2. 数据集和数学模型
  3. 回归参数估计
  4. 回归方程的显著性检验
  5. 残差分析和异常点检测
  6. 模型预测

1. 一元线性回归介绍

回归分析(Regression Analysis)是用来确定2个或2个以上变量间关系的一种统计分析方法。如果回归分析中,只包括一个自变量X和一个因变量Y时,且它们的关系是线性的,那么这种回归分析称为一元线性回归分析。

回归分析属于统计学的基本模型,涉及统计学基础,就会有一大堆的名词和知识点需要介绍。

在回归分析中,变量有2类:因变量 和 自变量。因变量通常是指实际问题中所关心的指标,用Y表示。而自变量是影响因变量取值的一个变量,用X表示,如果有多个自变量则表示为X1, X2, …, Xn。

回归分析研究的主要步骤:

  1. 确定因变量Y 与 自变量X1, X2, …, Xn 之间的定量关系表达式,即回归方程。
  2. 对回归方程的置信度检查。
  3. 判断自变量Xn(n=1,2,…,m)对因变量的影响。
  4. 利用回归方程进行预测。

本文会根据回归分析的的主要步骤,进行结构梳理,介绍一元线性回归模型的使用方法。

reg

2. 数据集和数学模型

先让我们通过一个例子开始吧,用一组简单的数据来说明一元线性回归分析的数学模型的原理和公式。找出下面数据集中Y与X的定量关系。

数据集为2016年3月1日,白天开盘的交易数据,为锌的2个期货合约的分钟线的价格数据。数据集包括有3列,索引列为时间,zn1.Close为ZN1604合约的1分钟线的报价数据,zn2.Close为ZN1605合约的1分钟线的报价数据。

数据集如下:


                    zn1.Close zn2.Close
2016-03-01 09:01:00     14075     14145
2016-03-01 09:02:00     14095     14160
2016-03-01 09:03:00     14095     14160
2016-03-01 09:04:00     14095     14165
2016-03-01 09:05:00     14120     14190
2016-03-01 09:06:00     14115     14180
2016-03-01 09:07:00     14110     14170
2016-03-01 09:08:00     14110     14175
2016-03-01 09:09:00     14105     14170
2016-03-01 09:10:00     14105     14170
2016-03-01 09:11:00     14120     14180
2016-03-01 09:12:00     14105     14170
2016-03-01 09:13:00     14105     14170
2016-03-01 09:14:00     14110     14175
2016-03-01 09:15:00     14105     14175
2016-03-01 09:16:00     14120     14185
2016-03-01 09:17:00     14125     14190
2016-03-01 09:18:00     14115     14185
2016-03-01 09:19:00     14135     14195
2016-03-01 09:20:00     14125     14190
2016-03-01 09:21:00     14135     14205
2016-03-01 09:22:00     14140     14210
2016-03-01 09:23:00     14140     14200
2016-03-01 09:24:00     14135     14205
2016-03-01 09:25:00     14140     14205
2016-03-01 09:26:00     14135     14205
2016-03-01 09:27:00     14130     14205

我们以zn1.Close列的价格为X,zn2.Close列的价格为Y,那么试试找到自变量X和因变量Y的关系的表达式。

为了直观起见,我们可以先画出一张散点图,以X为横坐标,Y为纵坐标,每个点对应一个X和一个Y。


# 数据集已存在df变量中
> head(df)
                    zn1.Close zn2.Close
2016-03-01 09:01:00     14075     14145
2016-03-01 09:02:00     14095     14160
2016-03-01 09:03:00     14095     14160
2016-03-01 09:04:00     14095     14165
2016-03-01 09:05:00     14120     14190
2016-03-01 09:06:00     14115     14180

# 分别给x,y赋值
> x<-as.numeric(df[,1])
> y<-as.numeric(df[,2])

# 画图
> plot(y~x+1)

01

从散点图上发现 X和Y 的排列基本是在一条直线附近,那么我们可以假设X和Y的关系是线性,可以用公式表式为。

Y = a + b * X + c
  • Y,为因变量
  • X,为自变量
  • a,为截距
  • b,为自变量系数
  • a+b*X, 表示Y随X的变化而线性变化的部分
  • c, 为残差或随机误差,是其他一切不确定因素影响的总和,其值不可观测。假定c是符合均值为0方差为σ^2的正态分布 ,记作c~N(0,σ^2)

对于上面的公式,称函数f(X) = a + b * X 为一元线性回归函数,a为回归常数,b为回归系数,统称回归参数。X 为回归自变量或回归因子,Y 为回归因变量或响应变量。如果(X1,Y1),(X2,Y2)…(Xn,Yn)是(X,Y)的一组观测值,则一元线性回归模型可表示为


Yi = a + b * X + ci,     i= 1,2,...n

其中E(ci)=0, var(ci)=σ^2, i=1,2,...n

通过对一元线性回归模型的数学定义,接下来让我们利用数据集做回归模型的参数估计。

3. 回归参数估计

对于上面的公式,回归参数a,b是我们不知道的,我们需要用参数估计的方法来计算出a,b的值,而从得到数据集的X和Y的定量关系。我们的目标是要计算出一条直线,使直接线上每个点的Y值和实际数据的Y值之差的平方和最小,即(Y1实际-Y1预测)^2+(Y2实际-Y2预测)^2+ …… +(Yn实际-Yn预测)^2 的值最小。参数估计时,我们只考虑Y随X的线性变化的部分,而残差c是不可观测的,参数估计法并不需要考虑残差,对于残差的分析在后文中介绍。

令公式变形为a和b的函数Q(a,b), 即 (Y实际-Y测试)的平方和,变成到(Y实际 – (a+b*X))的平方和。

reg2

公式一 回归参数变形公式

通过最小二乘估计推导出a和b的求解公式,详细的推导过程请参考文章一元线性回归的细节

reg3

公式二 回归参数计算公式

其中 x和y的均值,计算方法如下

reg4

公式三 均值计算公式

有了这个公式,我们就可以求出a和b两个的回归参数的解了。

接下来,我们用R语言来实现对上面数据的回归模型的参数估计,R语言中可以用lm()函数来实现一元线性回归的建模过程。


# 建立线性回归模型
> lm.ab<-lm(y ~ 1+x)

# 打印参数估计的结果
> lm.ab

Call:
lm(formula = y ~ 1 + x)

Coefficients:
(Intercept)            x  
   -349.493        1.029  

如果你想动手来计算也可以自己实现公式。


# x均值
> Xm<-mean(x);Xm 
[1] 14034.82

# y均值
> Ym<-mean(y);Ym
[1] 14096.76

# 计算回归系数
> b <- sum((x-Xm)*(y-Ym)) / sum((x-Xm)^2) ;b
[1] 1.029315

# 计算回归常数
> a <- Ym - b * Xm;a
[1] -349.493

回归参数a和b的计算结果,与lm()函数的计算结果是相同的。有了a和b的值,我们就可以画出这条近似的直接线。

计算公式为:

Y= a + b * X = -349.493 + 1.029315 * X 

画出回归线。


> plot(y~x+1)
> abline(lm.ab)

02

这条直线是我们用数据拟合出来的,是一个近似的值。我们看到有些点在线上,有些点不在线上。那么要评价这条回归线拟合的好坏,我们就需要对回归模型进行显著性检验。

4. 回归方程的显著性检验

从回归参数的公式二可知,在计算过程中并不一定要知道Y和X是否有线性相关的关系。如果不存相关关系,那么回归方程就没有任何意义了,如果Y和X是有相关关系的,即Y会随着X的变化而线性变化,这个时候一元线性回归方程才有意义。所以,我们需要用假设检验的方法,来验证相关性的有效性。

通常会采用三种显著性检验的方法。

  • T检验法:T检验是检验模型某个自变量Xi对于Y的显著性,通常用P-value判断显著性,小于0.01更小时说明这个自变量Xi与Y相关关系显著。
  • F检验法:F检验用于对所有的自变量X在整体上看对于Y的线性显著性,也是用P-value判断显著性,小于0.01更小时说明整体上自变量与Y相关关系显著。
  • R^2(R平方)相关系统检验法:用来判断回归方程的拟合程度,R^2的取值在0,1之间,越接近1说明拟合程度越好。

在R语言中,上面列出的三种检验的方法都已被实现,我们只需要把结果解读。上文中,我们已经通过lm()函数构建一元线性回归模型,然后可以summary()函数来提取模型的计算结果。


> summary(lm.ab)      # 计算结果

Call:
lm(formula = y ~ 1 + x)

Residuals:
     Min       1Q   Median       3Q      Max 
-11.9385  -2.2317  -0.1797   3.3546  10.2766 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -3.495e+02  7.173e+01  -4.872 2.09e-06 ***
x            1.029e+00  5.111e-03 201.390  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 4.232 on 223 degrees of freedom
Multiple R-squared:  0.9945,	Adjusted R-squared:  0.9945 
F-statistic: 4.056e+04 on 1 and 223 DF,  p-value: < 2.2e-16

模型解读:

  • Call,列出了回归模型的公式。
  • Residuals,列出了残差的最小值点,1/4分位点,中位数点,3/4分位点,最大值点。
  • Coefficients,表示参数估计的计算结果。
  • Estimate,为参数估计列。Intercept行表示常数参数a的估计值 ,x行表示自变量x的参数b的估计值。
  • Std. Error,为参数的标准差,sd(a), sd(b)
  • t value,为t值,为T检验的值
  • Pr(>|t|) ,表示P-value值,用于T检验判定,匹配显著性标记
  • 显著性标记,***为非常显著,**为高度显著, **为显著,·为不太显著,没有记号为不显著。
  • Residual standard error,表示残差的标准差,自由度为n-2。
  • Multiple R-squared,为相关系数R^2的检验,越接近1则越显著。
  • Adjusted R-squared,为相关系数的修正系数,解决多元回归自变量越多,判定系数R^2越大的问题。
  • F-statistic,表示F统计量,自由度为(1,n-2),p-value:用于F检验判定,匹配显著性标记。

通过查看模型的结果数据,我们可以发现通过T检验的截距和自变量x都是非常显著,通过F检验判断出整个模型的自变量是非常显著,同时R^2的相关系数检验可以判断自变量和因变量是高度相关的。

最后,我们通过的回归参数的检验与回归方程的检验,得到最后一元线性回归方程为:

Y = -349.493 + 1.029315 * X 

5. 残差分析和异常点检测

在得到的回归模型进行显著性检验后,还要在做残差分析(预测值和实际值之间的差),检验模型的正确性,残差必须服从正态分布N(0,σ^2)。

我们可以自己计算数据残差,并进行正态分布检验。


# 残差
> y.res<-residuals(lm.ab)

# 打印前6条数据
> head(y.res)
        1         2         3         4         5         6 
6.8888680 1.3025744 1.3025744 6.3025744 5.5697074 0.7162808 

# 正态分布检验
> shapiro.test(y.res)

	Shapiro-Wilk normality test

data:  y.res
W = 0.98987, p-value = 0.1164

# 画出残差散点图
> plot(y.res)

09

对残差进行Shapiro-Wilk正态分布检验,W接近1,p-value>0.05,证明数据集符合正态分布!关于正态分布的介绍,请参考文章常用连续型分布介绍及R语言实现

同时,我们也可以用R语言的工具生成4种用于模型诊断的图形,简化自己写代码计算的操作。


# 画图,回车展示下一张
> plot(lm.ab)    
Hit  to see next plot:   # 残差拟合图
Hit  to see next plot:   # 残差QQ图
Hit  to see next plot:   # 标准化的残差对拟合值 
Hit  to see next plot:   # 标准化残差对杠杆值

04

图1,残差和拟合值对比图

对残差和拟合值作图,横坐标是拟合值,纵坐标是残差。残差和拟合值之间,数据点均匀分布在y=0两侧,呈现出随机的分布,红色线呈现出一条平稳的曲线并没有明显的形状特征,说明残差数据表现非常好。

05

图2,残差QQ图

残差QQ图,用来描述残差是否符合正态分布。图中的数据点按对角直线排列,趋于一条直线,并被对角直接穿过,直观上符合正态分布。对于近似服从正态分布的标准化残差,应该有 95% 的样本点落在 [-2,2] 区间内。

06

图3,标准化残差平方根和拟合值对比图

对标准化残差平方根和拟合值作图,横坐标是拟合值,纵坐标是标准化后的残差平方根。与残差和拟合值对比图(图1)的判断方法类似,数据随机分布,红色线呈现出一条平稳的曲线,无明显的形状特征。

07

图4,标准残差和杠杆值对比图

对标准化残差和杠杆值作图,虚线表示的cooks距离等高线,通常用Cook距离度量的回归影响点。本图中没有出现红色的等高线,则说明数据中没有特别影响回归结果的异常点。

如果想把把4张图画在一起进行展示,可以改变画布布局。


> par(mfrow=c(2,2))
> plot(lm.ab)

08

看到上面4幅中,每幅图上都有一些点被特别的标记出来了,这些点是可能存在的异常值点,如果要对模型进行优化,我们可以从这些来入手。但终于本次残差分析的结果已经很好了,所以对于异常点的优化,可能并不能明显的提升模型的效果。

从图中发现,索引编号为27和192的2个点在多幅图中出现。我们假设这2个点为异常点,从数据中去掉这2个点,再进行显著性检验和残差分析。


# 查看27和192
> df[c(27,192),]
                    zn1.Close zn2.Close
2016-03-01 09:27:00     14130     14205
2016-03-01 14:27:00     14035     14085

# 新建数据集,去掉27和192
> df2<-df[-c(27,192),]

回归建模和显著性检验。


> x2<-as.numeric(df2[,1])
> y2<-as.numeric(df2[,2])
> lm.ab2<-lm(y2 ~ 1+x2)
> summary(lm.ab2)

Call:
lm(formula = y2 ~ 1 + x2)

Residuals:
    Min      1Q  Median      3Q     Max 
-9.0356 -2.1542 -0.2727  3.3336  9.5879 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -3.293e+02  7.024e+01  -4.688 4.83e-06 ***
x2           1.028e+00  5.004e-03 205.391  < 2e-16 ***
---
Signif. codes:  
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 4.117 on 221 degrees of freedom
Multiple R-squared:  0.9948,	Adjusted R-squared:  0.9948 
F-statistic: 4.219e+04 on 1 and 221 DF,  p-value: < 2.2e-16

对比这次的显著性检验结果和之前结果,T检验,F检验 和 R^2检验,并没有明显的效果提升,结果和我预想的是一样的。所以,通过残差分析和异常点分析,我认为模型是有效的。

6. 模型预测

最后,我们获得了一元线性回归方程的公式,就可以对数据进行预测了。比如,对给定X=x0时,计算出y0=a+b*x0的值,并计算出置信度为1-α的预测区间。

当X=x0,Y=y0时,置信度为1-α的预测区间为

reg5

reg6

我们可以用R语言的predict()函数来计算预测值y0,和相应的预测区间。程序算法如下。


> new<-data.frame(x=14040)
> lm.pred<-predict(lm.sol,new,interval="prediction",level=0.95)

# 预测结果
> lm.pred
       fit      lwr      upr
1 14102.09 14093.73 14110.44

当x0=14040时,在预测区间为0.95的概率时,y0的值为 14102,预测区间为[14093.73,14110.44]。

我们通过图形来表示。


> plot(y~x+1)
> abline(lm.ab,col='red')
> points(rep(newX$x,3),y=lm.pred,pch=19,col=c('red','blue','green'))

03

其中,红色点为y0的值,蓝色点为预测区间最小值,绿色点为预测区间最大值。

对于统计模型中最核心部分就在结果解读,本文介绍了一元回归模型的基本的建模过程和模型的详细解读方法。在我们掌握了这种方法以后,就可以更容易地理解和学习 多元回归,非线性回归 等更多的模型,并把这些模型应用到实际的工作中了。下一篇文章R语言解读多元线性回归模型,请继续阅读。

转载请注明出处:
http://blog.fens.me/r-linear-regression/

打赏作者

R语言中文分词包jiebaR

R的极客理想系列文章,涵盖了R的思想,使用,工具,创新等的一系列要点,以我个人的学习和体验去诠释R的强大。

R语言作为统计学一门语言,一直在小众领域闪耀着光芒。直到大数据的爆发,R语言变成了一门炙手可热的数据分析的利器。随着越来越多的工程背景的人的加入,R语言的社区在迅速扩大成长。现在已不仅仅是统计领域,教育,银行,电商,互联网….都在使用R语言。

要成为有理想的极客,我们不能停留在语法上,要掌握牢固的数学,概率,统计知识,同时还要有创新精神,把R语言发挥到各个领域。让我们一起动起来吧,开始R的极客理想。

关于作者:

  • 张丹(Conan), 程序员Java,R,PHP,Javascript
  • weibo:@Conan_Z
  • blog: http://blog.fens.me
  • email: bsspirit@gmail.com

转载请注明出处:
http://blog.fens.me/r-word-jiebar/

jiebaR

前言

本文挖掘是数据挖掘中一个非常重要的部分,有非常广阔的使用场景,比如我们可以对新闻事件进行分析,了解国家大事;也可以对微博信息进行分析,通过社交舆情看看大家的关注点。通过文本挖掘找到文章中的隐藏信息,对文章的结构进行分析,判断是不是同一个作者写文章;同时可以对邮件分析,结合bayes算法判断哪些是垃圾邮件,哪些是有用的邮件。

本文挖掘的第一步,就是要进行分词,分词将直接影响文本挖掘的效果。R语言在分词方面有很好的支持,接下来就给大家介绍一个不错的R语言中文分词包“结巴分词”(jiebaR)。

目录

  1. jiebaR包介绍
  2. 5分钟上手
  3. 分词引擎
  4. 配置词典
  5. 停止词过滤
  6. 关键词提取

1. jiebaR包介绍

结巴分词(jiebaR),是一款高效的R语言中文分词包,底层使用的是C++,通过Rcpp进行调用很高效。结巴分词基于MIT协议,就是免费和开源的,感谢国人作者的给力支持,让R的可以方便的处理中文文本。

官方Github的地址:https://github.com/qinwf/jiebaR

本文所使用的系统环境

  • Win10 64bit
  • R: 3.2.3 x86_64-w64-mingw32/x64 b4bit

jiebaR包是在CRAN发布的标准库,安装起来非常简单,2条命令就可以了。


~ R
> install.packages("jiebaR")
> library("jiebaR")

如果想要安装开发版本,可以使用devtools来进行安装,devtools的介绍请参考文章:在巨人的肩膀前行 催化R包开发


> library(devtools)
> install_github("qinwf/jiebaRD")
> install_github("qinwf/jiebaR")
> library("jiebaR")

开发版本安装,官方建议使用Linux系统 gcc >= 4.6 编译,Windows需要安装 Rtools。

2. 5分钟上手

5分钟上手,直接看第一个例子吧,对一段文字进行分词。


> wk = worker()

> wk["我是《R的极客理想》图书作者"]
[1] "我是" "R"    "的"   "极客" "理想" "图书" "作者"

> wk["我是R语言的深度用户"]
[1] "我"   "是"   "R"    "语言" "的"   "深度" "用户"

很简单地,2行代码,就完成了中文分词。

jiebaR提供了3种分词语句的写法,例子上面的用[]符号的语法,还可以使用<=符合语法,或者使用segment()函数。虽然形式不同,但是分词效果是一样的。 使用<=符号的语法,如下


> wk<='另一种符合的语法'
[1] "另"   "一种" "符合" "的"   "语法"

使用segment()函数的语法,如下


> segment( "segment()函数语句的写法" , wk )
[1] "segment" "函数"    "语句"    "的"      "写法" 

如果你觉得很神奇,想了解如何自定义操作符的,可以检查项目的源代码quick.R文件


# <= 符号定义
`<=.qseg`<-function(qseg, code){
  if(!exists("quick_worker",envir = .GlobalEnv ,inherits = F) || 
       .GlobalEnv$quick_worker$PrivateVarible$timestamp != TIMESTAMP){
    
    if(exists("qseg",envir = .GlobalEnv,inherits = FALSE ) ) 
      rm("qseg",envir = .GlobalEnv)
    
    modelpath  = file.path(find.package("jiebaR"),"model","model.rda")
    quickparam = readRDS(modelpath)
    
    if(quickparam$dict == "AUTO") quickparam$dict = DICTPATH
    if(quickparam$hmm == "AUTO") quickparam$hmm = HMMPATH
    if(quickparam$user == "AUTO") quickparam$user = USERPATH
    if(quickparam$stop_word == "AUTO") quickparam$stop_word = STOPPATH
    if(quickparam$idf == "AUTO") quickparam$idf = IDFPATH
    
    createquickworker(quickparam)
    setactive()
  } 

  //..代码省略
}

# [ 符号定义
`[.qseg`<- `<=.qseg`

我们也可以直接对文本文件进行分词,在当前目录新建一个文本文件idea.txt。


~ notepad idea.txt

R的极客理想系列文章,涵盖了R的思想,使用,工具,创新等的一系列要点,以我个人的学习和体验去诠释R的强大。

R语言作为统计学一门语言,一直在小众领域闪耀着光芒。直到大数据的爆发,R语言变成了一门炙手可热的数据分析的利器。随着越来越多的工程背景的人的加入,R语言的社区在迅速扩大成长。现在已不仅仅是统计领域,教育,银行,电商,互联网….都在使用R语言。

当然,我们运行分词程序,会在当前目录生成一个新的分词结果的文件。


> wk['./idea.txt']
[1] "./idea.segment.2016-07-20_23_25_34.txt"

打开文件idea.segment.2016-07-20_23_25_34.txt,整个本文以空格进行分词。


~ notepad idea.segment.2016-07-20_23_25_34.txt

R 的 极客 理想 系列 文章 涵盖 了 R 的 思想 使用 工具 创新 等 的 一系列 要点 以 我 个人 的 学习 和 体验 去 诠释 R 的 强大 R 语言 作为 统计学 一门 语言 一直 在 小众 领域 闪耀着 光芒 直到 大 数据 的 爆发 R 语言 变成 了 一门 炙手可热 的 数据分析 的 利器 随着 越来越 多 的 工程 背景 的 人 的 加入 R 语言 的 社区 在 迅速 扩大 成长 现在 已 不仅仅 是 统计 领域 教育 银行 电商 互联网 都 在 使用 R 语言

是不是很简单,5分钟实践就能完成分词的任务。

3. 分词引擎

在调用worker()函数时,我们实际是在加载jiebaR库的分词引擎。jiebaR库提供了7种分词引擎。

  • 混合模型(MixSegment):是四个分词引擎里面分词效果较好的类,结它合使用最大概率法和隐式马尔科夫模型。
  • 最大概率法(MPSegment) :负责根据Trie树构建有向无环图和进行动态规划算法,是分词算法的核心。
  • 隐式马尔科夫模型(HMMSegment):是根据基于人民日报等语料库构建的HMM模型来进行分词,主要算法思路是根据(B,E,M,S)四个状态来代表每个字的隐藏状态。 HMM模型由dict/hmm_model.utf8提供。分词算法即viterbi算法。
  • 索引模型(QuerySegment):先使用混合模型进行切词,再对于切出来的较长的词,枚举句子中所有可能成词的情况,找出词库里存在。
  • 标记模型(tag)
  • Simhash模型(simhash)
  • 关键词模型(keywods)

如果你不太关心引擎的事,那么直接用官方推荐的混合模型(默认选择)就行了。查看worker()函数的定义。


worker(type = "mix", dict = DICTPATH, hmm = HMMPATH, user = USERPATH,
  idf = IDFPATH, stop_word = STOPPATH, write = T, qmax = 20, topn = 5,
  encoding = "UTF-8", detect = T, symbol = F, lines = 1e+05,
  output = NULL, bylines = F, user_weight = "max")

参数列表:

  • type, 引擎类型
  • dict, 系统词典
  • hmm, HMM模型路径
  • user, 用户词典
  • idf, IDF词典
  • stop_word, 关键词用停止词库
  • write, 是否将文件分词结果写入文件,默认FALSE
  • qmax, 最大成词的字符数,默认20个字符
  • topn, 关键词数,默认5个
  • encoding, 输入文件的编码,默认UTF-8
  • detect, 是否编码检查,默认TRUE
  • symbol, 是否保留符号,默认FALSE
  • lines, 每次读取文件的最大行数,用于控制读取文件的长度。大文件则会分次读取。
  • output, 输出路径
  • bylines, 按行输出
  • user_weight, 用户权重

我们在调用worker()时,就加载了分词引擎,可以打印出来,查看分词的引擎的配置。


> wk = worker()
> wk
Worker Type:  Jieba Segment

Default Method  :  mix     # 混合模型
Detect Encoding :  TRUE    # 检查编码
Default Encoding:  UTF-8   # UTF-8
Keep Symbols    :  FALSE   # 不保留符号
Output Path     :          # 输出文件目录
Write File      :  TRUE    # 写文件
By Lines        :  FALSE   # 不行输出
Max Word Length :  20      # 最大单单词长度
Max Read Lines  :  1e+05   # 最大读入文件行数

Fixed Model Components:  

$dict                      # 系统词典
[1] "D:/tool/R-3.2.3/library/jiebaRD/dict/jieba.dict.utf8"

$user                      # 用户词典
[1] "D:/tool/R-3.2.3/library/jiebaRD/dict/user.dict.utf8"

$hmm                       # 隐式马尔科夫模型模型
[1] "D:/tool/R-3.2.3/library/jiebaRD/dict/hmm_model.utf8"

$stop_word                 # 停止词,无
NULL

$user_weight               # 用户词典权重
[1] "max"

$timestamp                 # 时间戳
[1] 1469027302

$default $detect $encoding $symbol $output $write $lines $bylines can be reset.

如果我们想改变分词引擎的配置项,可以在调用worker()创建分词引擎时,也可以通过wk$XX来进行设置。如果想了解wk是什么类型的对象,我们通过pryr包的otype的函数来检查wk对象的类型。关于pryr包的详细使用,请参考文章撬动R内核的高级工具包pryr


# 加载 pryr包
> library(pryr)
> otype(wk)  # 面向对象的类型检查
[1] "S3"

> class(wk)  # 查看class是属性
[1] "jiebar"  "segment" "jieba" 

4. 配置词典

对于分词的结果好坏的关键因素是词典,jiebaR默认有配置标准的词典。对于我们的使用来说,不同行业或不同的文字类型,最好用专门的分词词典。在jiebaR中通过show_dictpath()函数可以查看默认的标准词典,可以通过上一小节介绍的配置项,来指定我们自己的词典。日常对话的常用词典,比如搜狗输入法的词库。


# 查看默认的词库位置
> show_dictpath()
[1] "D:/tool/R-3.2.3/library/jiebaRD/dict"

# 查看目录
> dir(show_dictpath())
[1] "D:/tool/R-3.2.3/library/jiebaRD/dict"
 [1] "backup.rda"      "hmm_model.utf8"  "hmm_model.zip"  
 [4] "idf.utf8"        "idf.zip"         "jieba.dict.utf8"
 [7] "jieba.dict.zip"  "model.rda"       "README.md"      
[10] "stop_words.utf8" "user.dict.utf8" 

看到词典目录中,包括了多个文件。

  • jieba.dict.utf8, 系统词典文件,最大概率法,utf8编码的
  • hmm_model.utf8, 系统词典文件,隐式马尔科夫模型,utf8编码的
  • user.dict.utf8, 用户词典文件,utf8编码的
  • stop_words.utf8,停止词文件,utf8编码的
  • idf.utf8,IDF语料库,utf8编码的
  • jieba.dict.zip,jieba.dict.utf8的压缩包
  • hmm_model.zip,hmm_model.utf8的压缩包
  • idf.zip,idf.utf8的压缩包
  • backup.rda,无注释
  • model.rda,无注释
  • README.md,说明文件

打开系统词典文件jieba.dict.utf8,并打印前50行。


> scan(file="D:/tool/R-3.2.3/library/jiebaRD/dict/jieba.dict.utf8",
+           what=character(),nlines=50,sep='\n',
+           encoding='utf-8',fileEncoding='utf-8')
Read 50 items
 [1] "1号店 3 n"  "1號店 3 n"  "4S店 3 n"   "4s店 3 n"  
 [5] "AA制 3 n"   "AB型 3 n"   "AT&T 3 nz"  "A型 3 n"   
 [9] "A座 3 n"    "A股 3 n"    "A輪 3 n"    "A轮 3 n"   
[13] "BB机 3 n"   "BB機 3 n"   "BP机 3 n"   "BP機 3 n"  
[17] "B型 3 n"    "B座 3 n"    "B股 3 n"    "B超 3 n"   
[21] "B輪 3 n"    "B轮 3 n"    "C# 3 nz"    "C++ 3 nz"  
[25] "CALL机 3 n" "CALL機 3 n" "CD机 3 n"   "CD機 3 n"  
[29] "CD盒 3 n"   "C座 3 n"    "C盘 3 n"    "C盤 3 n"   
[33] "C語言 3 n"  "C语言 3 n"  "D座 3 n"    "D版 3 n"   
[37] "D盘 3 n"    "D盤 3 n"    "E化 3 n"    "E座 3 n"   
[41] "E盘 3 n"    "E盤 3 n"    "E通 3 n"    "F座 3 n"   
[45] "F盘 3 n"    "F盤 3 n"    "G盘 3 n"    "G盤 3 n"   
[49] "H盘 3 n"    "H盤 3 n"

我们发现系统词典每一行都有三列,并以空格分割,第一列为词项,第二列为词频,第三列为词性标记。

打开用户词典文件user.dict.utf8,并打印前50行。


> scan(file="D:/tool/R-3.2.3/library/jiebaRD/dict/user.dict.utf8",
+      what=character(),nlines=50,sep='\n',
+      encoding='utf-8',fileEncoding='utf-8')
Read 5 items
[1] "云计算"   "韩玉鉴赏" "蓝翔 nz"  "CEO"      "江大桥"  

用户词典第一行有二列,,第一列为词项,第二列为词性标记,没有词频的列。用户词典默认词频为系统词库中的最大词频。

jiebaR包关于词典词性标记,采用ictclas的标记方法。ICTCLAS 汉语词性标注集。

代码 名称 帮助记忆的诠释
Ag 形语素 形容词性语素。形容词代码为a,语素代码g前面置以A。
a 形容词 取英语形容词adjective的第1个字母。
ad 副形词 直接作状语的形容词。形容词代码a和副词代码d并在一起。
an 名形词 具有名词功能的形容词。形容词代码a和名词代码n并在一起。
b 区别词 取汉字"别"的声母。
c 连词 取英语连词conjunction的第1个字母。
Dg 副语素 副词性语素。副词代码为d,语素代码g前面置以D。
d 副词 取adverb的第2个字母,因其第1个字母已用于形容词。
e 叹词 取英语叹词exclamation的第1个字母。
f 方位词 取汉字"方"的声母。
g 语素 绝大多数语素都能作为合成词的"词根",取汉字"根"的声母。
h 前接成分 取英语head的第1个字母。
i 成语 取英语成语idiom的第1个字母。
j 简称略语 取汉字"简"的声母。
k 后接成分
l 习用语 习用语尚未成为成语,有点"临时性",取"临"的声母。
m 数词 取英语numeral的第3个字母,n,u已有他用。
Ng 名语素 名词性语素。名词代码为n,语素代码g前面置以N。
n 名词 取英语名词noun的第1个字母。
nr 人名 名词代码n和"人(ren)"的声母并在一起。
ns 地名 名词代码n和处所词代码s并在一起。
nt 机构团体 "团"的声母为t,名词代码n和t并在一起。
nz 其他专名 "专"的声母的第1个字母为z,名词代码n和z并在一起。
o 拟声词 取英语拟声词onomatopoeia的第1个字母。
p 介词 取英语介词prepositional的第1个字母。
q 量词 取英语quantity的第1个字母。
r 代词 取英语代词pronoun的第2个字母,因p已用于介词。
s 处所词 取英语space的第1个字母。
Tg 时语素 时间词性语素。时间词代码为t,在语素的代码g前面置以T。
t 时间词 取英语time的第1个字母。
u 助词 取英语助词auxiliary 的第2个字母,因a已用于形容词。
Vg 动语素 动词性语素。动词代码为v。在语素的代码g前面置以V。
v 动词 取英语动词verb的第一个字母。
vd 副动词 直接作状语的动词。动词和副词的代码并在一起。
vn 名动词 指具有名词功能的动词。动词和名词的代码并在一起。
w 标点符号
x 非语素字 非语素字只是一个符号,字母x通常用于代表未知数、符号。
y 语气词 取汉字"语"的声母。
z 状态词 取汉字"状"的声母的前一个字母。

下面我们自定义一个用户词典,来试试效果。编写词典文件,user.utf8。


~ notepad user.utf8

R语言
R的极客理想
大数据
数据

使用我们的自定义的用户词典,对刚才的文本再进行分词。


> wk = worker(user='user.utf8')
> wk['./idea.txt']
[1] "./idea.segment.2016-07-21_11_14_24.txt"

对比2次产生的分词结果,idea.segment.2016-07-20_23_25_34.txt 和 idea.segment.2016-07-21_11_14_24.txt。

jiebaR-cut

在实际使用中,jiebaR默认提供的用户词典只有5个单词,太简单了,肯定是不够用的。我们可以用搜狗词典,来丰富用户自己的词库。接下来,让我们配置搜狗词典。你需要安装一个搜狗输入法,具体的安装过程不再解释。

我安装的是搜狗五笔输入法,找到搜狗的安装目录,并找到词典文件。我的搜狗词典,在下面的安装位置。


C:\Program Files (x86)\SogouWBInput\2.1.0.1288\scd\17960.scel

把17960.scel文件复制到自己的项目目录里,用文本编辑器打开文件,发现是二进制的。那么我需要用工具进行转换,把二进制的词典转成我们可以使用的文本文件。jiebaR包的作者,同时开发了一个cidian项目,可以转换搜狗的词典,那么我们只需要安装cidian包即可。

安装cidian项目


> install.packages("devtools")
> install.packages("stringi")
> install.packages("pbapply")
> install.packages("Rcpp")
> install.packages("RcppProgress")
> library(devtools)
> install_github("qinwf/cidian")
> library(cidian)

转换二进制词典到文本文件。


# 转换
> decode_scel(scel = "./17960.scel",cpp = TRUE)
output file: ./17960.scel_2016-07-21_00_22_11.dict

# 查看生成的词典文件
> scan(file="./17960.scel_2016-07-21_00_22_11.dict",
+      what=character(),nlines=50,sep='\n',
+      encoding='utf-8',fileEncoding='utf-8')
Read 50 items
 [1] "阿坝州 n"         "阿百川 n"         "阿班 n"          
 [4] "阿宾 n"           "阿波菲斯 n"       "阿不都热希提 n"  
 [7] "阿不都西库尔 n"   "阿不力克木 n"     "阿尔姆格伦 n"    
[10] "阿尔沙文 n"       "阿肥星 n"         "阿菲正传 n"      
[13] "阿密特 n"         "阿穆 n"           "阿穆隆 n"        
[16] "阿帕鲁萨镇 n"     "阿披实 n"         "阿衰 n"          
[19] "阿霞 n"           "艾奥瓦 n"         "爱不疚 n"        
[22] "爱的错位 n"       "爱得得体 n"       "爱的火焰 n"      
[25] "爱的流刑地 n"     "爱得起 n"         "埃夫隆 n"        
[28] "爱搞网 n"         "爱国红心 n"       "爱呼 n"          
[31] "爱就宅一起 n"     "埃克希儿 n"       "爱没有错 n"      
[34] "埃蒙斯 n"         "爱奴新传 n"       "爱起点 n"        
[37] "爱情的牙齿 n"     "爱情海滨 n"       "爱情节 n"        
[40] "爱情美的样子 n"   "爱情无限谱 n"     "爱情占线 n"      
[43] "爱情转移 n"       "爱情左灯右行 n"   "爱上你是一个错 n"
[46] "矮哨兵 n"         "爱是妥协 n"       "爱似水仙 n"      
[49] "爱太痛 n"         "爱无界 n"    

接下来,直接把搜狗词典配置到我们的分词库中,就可以直接使用了。把搜狗词典文件改名,从17960.scel_2016-07-21_00_22_11.dict到user.dict.utf8,然后替换D:\tool\R-3.2.3\library\jiebaRD\dict目录下面的user.dict.utf8。这样默认的用户词典,就是搜狗词典了。很酷吧!

5. 停止词过滤

停止词就是分词过程中,我们不需要作为结果的词,像英文的语句中有很多的a,the,or,and等,中文语言中也有很多,比如 的,地,得,我,你,他。这些词因为使用频率过高,会大量出现在一段文本中,对于分词后的结果,在统计词频的时候会增加很多的噪音,所以我们通常都会将这些词进行过滤。

在jiebaR中,过滤停止词有2种方法,一种是通过配置stop_word文件,另一种是使用filter_segment()函数。

首先我们先来看,通过配置stop_word文件的方法。新建一个stop_word.txt文件。


~ notepad stop_word.txt

我
我是

加载分词引擎,并配置停止词过滤。


> wk = worker(stop_word='stop_word.txt')
> segment<-wk["我是《R的极客理想》图书作者"]
> segment
[1] "R"    "的"   "极客" "理想" "图书" "作者"

上面的文本,我们把"我是"通过停止词进行了过滤。如果还想过滤“作者”一词,可以动态的调用filter_segment()函数。


> filter<-c("作者")
> filter_segment(segment,filter)
[1] "R"    "的"   "极客" "理想" "图书"

6. 关键词提取

关键词提取是文本处理非常重要的一个环节,一个经典算法是TF-IDF算法。其中,TF(Term Frequency)代表词频,IDF(Inverse Document Frequency)表示逆文档频率。如果某个词在文章中多次出现,而且不是停止词,那么它很可能就反应了这段文章的特性,这就是我们要找的关键词。再通过IDF来算出每个词的权重,不常见的词出现的频率越高,则权重越大。计算TF-IDF的公式为:

TF-IDF = TF(词频) * 逆文档频率(IDF)

对文档中每个词计算TF-IDF的值,把结果从大到小排序,就得到了这篇文档的关键性排序列表。关于IF-IDF的解释,参考了文章TF-IDF与余弦相似性的应用(一):自动提取关键词

jiebaR包的关键词提取提取的实现,也是使用了TF-IDF的算法。在安装目录中的idf.utf8文件,为IDF的语料库。查看idf.utf8内容。


> scan(file="D:/tool/R-3.2.3/library/jiebaRD/dict/idf.utf8",
+      what=character(),nlines=50,sep='\n',
+      encoding='utf-8',fileEncoding='utf-8')
Read 50 items
 [1] "劳动防护 13.900677652"      "生化学 13.900677652"       
 [3] "奥萨贝尔 13.900677652"      "考察队员 13.900677652"     
 [5] "岗上 11.5027823792"         "倒车档 12.2912397395"      
 [7] "编译 9.21854642485"         "蝶泳 11.1926274509"        
 [9] "外委 11.8212361103"         "故作高深 11.9547675029"    
[11] "尉遂成 13.2075304714"       "心源性 11.1926274509"      
[13] "现役军人 10.642581114"      "杜勃留 13.2075304714"      
[15] "包天笑 13.900677652"        "贾政陪 13.2075304714"      
[17] "托尔湾 13.900677652"        "多瓦 12.5143832909"        
[19] "多瓣 13.900677652"          "巴斯特尔 11.598092559"     
[21] "刘皇帝 12.8020653633"       "亚历山德罗夫 13.2075304714"
[23] "社会公众 8.90346537821"     "五百份 12.8020653633"      
[25] "两点阈 12.5143832909"       "多瓶 13.900677652"         
[27] "冰天 12.2912397395"         "库布齐 11.598092559"       
[29] "龙川县 12.8020653633"       "银燕 11.9547675029"        
[31] "历史风貌 11.8212361103"     "信仰主义 13.2075304714"    
[33] "好色 10.0088573539"         "款款而行 12.5143832909"    
[35] "凳子 8.36728816325"         "二部 9.93038573842"        
[37] "卢巴 12.1089181827"         "五百五 13.2075304714"      
[39] "畅叙 11.598092559"          "吴栅子 13.2075304714"      
[41] "智力竞赛 13.900677652"      "库邦 13.2075304714"        
[43] "非正义 11.3357282945"       "编订 10.2897597393"        
[45] "悲号 12.8020653633"         "陈庄搭 13.2075304714"      
[47] "二郎 9.62401153296"         "电光石火 11.8212361103"    
[49] "抢球 11.9547675029"         "南澳大利亚 10.9562386728"  

idf.utf8文件每一行有2列,第一列是词项,第二列为权重。然后,我通过计算文档的词频(TF),与语料库的IDF值相乘,就可以得到TF-IDF值,从而提取文档的关键词。

比如,我们对下面的文本内容进行关键词的提取。


> wk = worker()
> segment<-wk["R的极客理想系列文章,涵盖了R的思想,使用,工具,创新等的一系列要点,以我个人的学习和体验去诠释R的强大。"]

# 计算词频
> freq(segment)
     char freq
1    创新    1
2      了    1
3    文章    1
4    强大    1
5       R    3
6    个人    1
7      的    5
8    诠释    1
9      和    1
10 一系列    1
11   使用    1
12     以    1
13     等    1
14   极客    1
15   理想    1
16   思想    1
17   涵盖    1
18   系列    1
19     去    1
20     我    1
21   工具    1
22   学习    1
23   体验    1
24   要点    1

# 取TF-IDF的前5的关键词
> keys = worker("keywords",topn=5)

# 计算关键词
> vector_keywords(segment,keys)
11.7392 8.97342 8.23425  8.2137 7.43298 
 "极客"  "诠释"  "要点"  "涵盖"  "体验" 

使用jiebaR包处理分词确实简单,几行的代码就能实现分词的各种算法操作。有了这个工具,我们就可以文档中,发现各种语言规则进行文本挖掘了。下篇文章让我们挖掘一下上市公司的公告吧,说不定能发现什么市场规则。

本文只是抛砖引玉地介绍了jiebaR包的使用方法,详细使用操作,请参考包作者的官方介绍。再次感谢jiebaR作者@qinwenfeng,为R语言在中文分词中提供了一套非常不错的工具包!

转载请注明出处:
http://blog.fens.me/r-word-jiebar/

打赏作者

2016天善智能交流会第22场: R语言为量化而生

跨界知识聚会系列文章,“知识是用来分享和传承的”,各种会议、论坛、沙龙都是分享知识的绝佳场所。我也有幸作为演讲嘉宾参加了一些国内的大型会议,向大家展示我所做的一些成果。从听众到演讲感觉是不一样的,把知识分享出来,你才能收获更多。

关于作者

  • 张丹(Conan), 程序员Java,R,Nodejs
  • weibo:@Conan_Z
  • blog: http://blog.fens.me
  • email: bsspirit@gmail.com

转载请注明出处:
http://blog.fens.me/meeting-hellobi-20160701/

meeting-hellobi

前言

感谢天善智能社区的邀请,有幸参加每周一期的跟数据有关的行业、工具、技术的交流盛宴,活动的口号是“Friday BI Fly 周五BI飞起来”。

目录

  1. 我的分享主题:R语言为量化而生
  2. 会议体验
  3. 自由讨论

1. 我的分享主题:R语言为量化而生

本次分享的主题 R语言为量化而生,主要内容来自我的一篇博客文章:R语言为量化而生。希望能够解释清楚,在量化投资中为什么要用R语言。从程序员的角度看,C++,Java,Python, C#都是可行方案;从数据人员的角度看,Excel, SAS, Matlab更是不错的。那么为什么是R语言呢,R语言的优势在哪里体现?

这类的问题,总是会被问到。那么答案,就在于你对量化这件事情的了解,和对各种编程语言的理解。最近3年,互联网在量化领域的大发展,以Quantopian为代表的在线策略研发平台,用Python做为核心语言,国内同样支持Python的平台也有 优矿聚宽米筐。这些平台主是面向程序员群体的平台,希望通过挖掘草根明星,来推动量化的发展。传统的量化交易软件,像文华MC, TB, TS 都有自己一套的脚本化的编程语言。有实力的专业团队,通常会自成体系的独立开发一套自己的系统。如果面向更广泛的人群,最常用的方法就是Wind导数据,Excel中拉个表出来。

所以,其实用什么语言不重要,关键是怎么理解做量化这件事情。那么R语言的天生优势就是数学计算,数据处理,免费开源,大量支持库。试试吧,你一定会喜欢的。

2. 会议体验

本次分享受天善智能社区的邀请,我真的非常高兴。天善智能是新一代的商业智能和大数据的垂直社区,聚集了大量的数据分析从业人员。活动介绍,https://ask.hellobi.com/blog/tianshansoft/4229。 本次活动同时有30个微信群进行直播,参加的人员,至少有2000人以上。可以天善智能社区,在行业的影响力是非常大的。

发个截图,体会一下微信同步直播的震撼吧!

wx

本此的分享基于微信的直播,我也第一次体验,要用纯文字的方式来进行介绍。想把一个事情说清楚,又增加了不少的难度。由于不能分享屏幕,代码部分会通过图片截屏。

本次活动的总结,https://ask.hellobi.com/blog/tianshansoft/4271,感谢天善社区的工作人员进行整理。

远程分享,就是没能与大家合照,有点遗憾!!贴张自己的照片吧。

01

3. 自由讨论

分享后,很多朋友都对于R语言都是非常的好奇,提了很多的问题,用户的参与性非常强。下列直接贴出用户的问题和我的回复。

1、替新手问一个,请教一下,R语言的数据分析应该从哪方面入手练习啊?因为目前工作上不是用R的,看完书之后想具体去试一下。

张丹: R其实上手很快,找一本书,认真操作练习一遍就上手了。

2、玉琴:不建议用for loop的原因是考虑到性能问题吗

张丹:for loop是调用的R的循环库,apply是调用C的循环库,性能差距还是很大的

3、来自20群的提问:提个问题,微软对R的收购会对R语言的发展产生什么影响?

张丹:我觉得这是正向发展的,是好事情。大公司看到了R的潜力!

4、尚林栋:R语言金融建模的具体步骤能说一下吗

金融建模的具体步骤,你可以参考这篇文章,http://blog.fens.me/finance-stock-ma/

5、刘嘉丰Alan:丹哥,现在有很多量化平台,提供打包好的函数,在线回测,和自己造轮子拿R语言相比,您觉得各有什么优势呢?

张丹:R的优势就是在数学计算,数据处理上。行业标准还没有统一,所以不一定在线平台的轮子就一定好用。但另外,我们从开发或使用的角度,更多的用到的R包,都是RStudio公司的产品,我觉得是RStudio在推动R的整个的进化过程。

6、我也觉得r语言不错,但经常想不到商业场景,到现在,我只是用它统计考勤,各种绩效kpi,每月算一次奖金,已经这样过去2年了,r语言路在何方哪?

张丹: 你所说的统计,只能说简单计数。比如,你要预测下个月的考勤情况,从而设计预算方案。你可能就需要做个回归分析,这时R就能给你很大的帮助了。生活和工作中,随处都是数据分析的场景。

7、Allen:r在拟合上感觉比python用起来更爽一些,其返回的结果较多

张丹:那么R和python比,R更面向数据,特别是对于没有编程基础的人。PYTHON,还是程序语言,还要了解程序结构,程序架构,代码量不会少。

有IT背景程序员,可能更倾向于PYTHON;如果没有IT背景,R更容易上手。

8、越中女儿:请教一个问题:quantmod对美股的实时接口很好用,对A股不支持,且A股基本面数据才更新到2013.09,请问有好用的ETL包么,类似于python的tushare那样对A股友好的,各种etl啊清洗的脏活累活感觉python更好啊,R就是安安静静做做统计,玩玩图形。

张丹: quantmod使用的是yahoo等国外的数据源,这些数据源本身没有A股数据,如果需要A股数据,用tushare还是不错的。 R特有的data.frame,matrix 等类型和操作方法,在python也需要单独去实现。

9、柠檬味的香草:最近想研究一些互联网文本数据与指数或各股走势的关系,但是在使用R语言处理文本数据不是很方便,丹哥可有一些强大的library推荐,对于非结构,文本数据的处理。

张丹:“尽量使用向量计算或矩阵计算的计算方法”,可以这样理解,对于一个二维结构,for需要2次,0(N^2)的时间复杂度。如果我们把数据,直接就按矩阵存储, 你让矩阵里的每个点都加1, 只需要算一次。Hadley提供的包,源代码我都看过,写很棒,也很实用。

r在拟合上感觉比python用起来更爽一些,其返回的结果较多

其实R有很多的第三方的包,已经有了大量的算法包,而其他语言相对较少。只是我们平时接触的不多,所以觉得用不到。R有大量的统计包,你可以从官方网站找到,输出的结果,大部分也都是统计的结果。

R所支持的行业领域,非常广泛。而工程的语言,不会做细粒度的区分,只是通用的解决方法。

10、郑州—金融数据:python有pandas.DataFrame,pandas应该是第三方的数据库结构吧?R的data.frame是内置的。

张丹:pandas.DataFrame,在底层处理,还需要对原PYTHON的数据结构做映射。当然他可以解决的很好,但你看到的内存结构,可能并不是真正的内存结构。

R内置数据类型,就可以理解是内存结构。不需要再考虑转换了。找一个自己熟悉的语言,大多数的功能,每种语言都是能实现。只有很细的领域,才会进一步区分。

11、RHaoop采用分布式并行计算,那请问如何解决需要嵌套循环的算法。

张丹:对于基于hadoop大数据的MR计算,建议做数学变成,通过数学的角度处理。我写过2个例子,一个是pagerank, 一个是itemcf。

12、@柠檬味的香草:想听听丹哥对传统数据挖掘转量化投资的建议。比如前景?竞争力?

张丹:量化投资,其实是IT人都想转的行业。你写的代码,不是通过工资来赚钱,而直接通过交易赚钱,代码的效用是最大化的。但这个行业竞争很大,聪明人都在这里,要么你的技术牛,要么你了解市场,要么的算法是独特的,不然也很难。

JhT: 做量化交易和策略的都是高智商的

越中女儿:我觉得量化对金融市场的理解比对技术本身更重要,R的需求应该会很快凸显出来。因为数据基础都有了,后面就是差会分析的人了。通常懂数据分析的程序员,比纯程序员待遇高。

13、老师,有好的spark或者hadoop入门的书吗,计算机能力弱和java不懂啊

张丹:hadoop有很多书了,我当初看的是 权威指南。spark的书不了解,我的是网上文档。

14、@Mia.W 学RHadoop需要对Hadoop或Mapreduce了解到什么程度,需要从头学hadoop或java吗

张丹:hadoop的MR的原理要了解,找到懂JAVA的同事,帮你把环境搭好。

15、@JhT 我是刚进来的,R的优势是什么?

张丹:R是免费开源的,CRAN上有8000多个包,遍布各行各业。R语言的3个特性,数学计算,数据建模,可视化。

16、@郑州—金融数据个人感觉商业上matlab比R和python支持度都要好,不管是分析,统计,挖掘还是量化方便,收费的毕竟是收费的

张丹:有商业推动,当然要比免费的好了。不过,像SAS和Matlab也在打通和R的接口,毕竟由全球第三方贡献包,要比一家公司提供的包要多很多的。

17、@越中女儿 有用R做过实盘风控么

张丹:有做,其实不太复杂。你把需要的实时数据,都同步存到redis中,用R在秒级调reids取数据,计算完成再写回去。

18、@Jason.k计算机8g内存,数据虽然行数不多,但是很多列,所以数据csv格式大小会高达几个G,这个规模数据量,内存应该是不够的。

张丹:R的机制,会把数据一次性加载到内存中。就算能读到内存,每次计算时,也会有中间变量,所以你的基础内存是不够的。而且对于win性能会更差。

最后,再次感谢 天善社区的小伙伴们的努力,谢谢大家!

转载请注明出处:
http://blog.fens.me/meeting-hellobi-20160701/

打赏作者

R语言为量化而生

R的极客理想系列文章,涵盖了R的思想,使用,工具,创新等的一系列要点,以我个人的学习和体验去诠释R的强大。

R语言作为统计学一门语言,一直在小众领域闪耀着光芒。直到大数据的爆发,R语言变成了一门炙手可热的数据分析的利器。随着越来越多的工程背景的人的加入,R语言的社区在迅速扩大成长。现在已不仅仅是统计领域,教育,银行,电商,互联网….都在使用R语言。

要成为有理想的极客,我们不能停留在语法上,要掌握牢固的数学,概率,统计知识,同时还要有创新精神,把R语言发挥到各个领域。让我们一起动起来吧,开始R的极客理想。

关于作者:

  • 张丹(Conan), 程序员Java,R,Nodejs
  • weibo:@Conan_Z
  • blog: http://blog.fens.me
  • email: bsspirit@gmail.com

转载请注明出处:
http://blog.fens.me/r-finance/

r-finance

前言

做数据分析的朋友,一定听说过R语言。R语言是一门统计语言,在数据分析领域优势是非常明显的。

本文以 “R语言,为量化而生”为题,说明R语言真的很适合做金融做量化策略。金融本身是玩数据行业,R的最大的优势就是数据分析,所以用R来做量化投资的策略,真是很配,不仅顺手而且方便,用了你就会知识。

本文将由3个方面来介绍,R语言做量化是多么的适合。

目录

  1. 为什么是R语言?
  2. R语言的数据处理和时间序列
  3. R语言和金融模型

1. 为什么是R语言?

那么为什么是R语言,而不是其他的语言? 先简单介绍一下,我们的个人经历。

我是一个程序员,从2004年开始接触Java写了10多年的Java程序,期间还尝试过多种编程语言,VB、PHP、Python、SAS、R、Nodejs,最后把自己锁定在R,Nodejs和Java。谈不上对每一种语言都有很深的理解,但是每种语言的特点还是有点心得。

之所以选择R,Nodejs和Java这3种语言,有一部分情怀,更多的是理性。从技术发展来看,编程开发变得越来越简单,10年前用JavaEE做一个简单的web项目至少要2人月,现在用Nodejs新人边学边搞只需10人天。而且随着业务的多样化,单一的技术已经不足以支撑业务的发展,业务在从传统的软件开发向互联网和数据产品的方向在进化。根据不同语言的特点,每种都将在开发中占据一席之地,而很难在出现一种语言统一天下的情况。

R语言将在数据分析领域发挥着重要的作用。R语言的3个特性,数学计算、数据建模和数据可视化。R语言封装了多种基础学科的计算函数,我们在R语言编程的过程中只需要调用这些计算函数,就可以构建出面向不同领域、不同业务的、复杂的数学模型。

另外,R的知识体系结构是复杂的,要想学好R,就必须把多学科的知识综合运用,而最大的难点不在于R语言本身,在于使用者的知识基础和综合运用的能力。

r-basic

图中我将R语言知识体系结构分为3个部分:IT技术 + 业务知识 + 基础学科。

  • IT技术:是数据大发展时代必备的技术之一,R语言就是我们应该要掌握的一门技术。
  • 业务知识:是市场经验和法则,不管你在什么公司,你都了解业务是什么,产品是什么,用户是谁,公司的价值在哪里!
  • 基础学科:是我们在学校里学到的理论知识,虽然当初学的时候并不理解,工作中如果你还能掌握并实际运用,那么这将是你最有价值的竞争力。

关于R的知识体系,可以参考文章,R语言知识体系概览

对于金融量化投资来说,刚好是一个交叉学科,你需要懂IT技术,熟悉金融市场的规则,有数学建模的能力。R语言,正好可以帮我们来解决这样的问题,所以“R语言,为量化而生”!

对于做过数据分析的人来说,大家都了解什么是最费时间的!!无疑就是数据处理的部分。

2. R语言的数据处理和时间序列

第二部分,我们来介绍一下R语言的数据类型和数据处理的一些方法。当然,本文并没有介绍如何入门R语言,新手入门请参考文章R的极客理想系列文章

2.1 基本数据类型

在R语言中,数据类型包括向量类型,字符串类型,数字类型,布尔类型,矩阵类型,数据框类型,list类型等,通常我们在使用R语言里做数据处理的时候,大部分都会以数据框(data.frame)类型为一个主要的数据内存类型来使用。

数据框(data.frame)类型是R语言内置的一种数据类型,我们可以简单地把它理解为,与关系型数据库中表的结构是类似的,是一种二维的数据结构。


# 新建一个数据框
> data.frame(A=1:6,B=LETTERS[1:6])
  A B
1 1 A
2 2 B
3 3 C
4 4 D
5 5 E
6 6 F

正是由于R语言内置了这样的数据类型,使我们从数据库读取数据或导入CSV格式的数据时,与R语言有了一个很好的映射关系,直接加载到R语言的内存中变成标准化数据格式。

然后,就可以基于标准化的数据格式,用R语言的功能函数来处理数据了。比如,对于做数据库开发的人员来说,他可以使用sqldf包,在R语言中通过SQL语句对数据进行数据变换。同时,也可以按着数据框(data.frame)的标准方法进行数据处理,通过约定的向量索引下标的方式来按行按列来读取数据,或使用功能函数处理数据。


# sqldf包的使用
> library(sqldf)
> sqldf('select * from iris limit 6')
  Sepal.Length Sepal.Width Petal.Length Petal.Width Species
1          5.1         3.5          1.4         0.2  setosa
2          4.9         3.0          1.4         0.2  setosa
3          4.7         3.2          1.3         0.2  setosa
4          4.6         3.1          1.5         0.2  setosa
5          5.0         3.6          1.4         0.2  setosa
6          5.4         3.9          1.7         0.4  setosa

# 向量索引
> iris[1:6,]
  Sepal.Length Sepal.Width Petal.Length Petal.Width Species
1          5.1         3.5          1.4         0.2  setosa
2          4.9         3.0          1.4         0.2  setosa
3          4.7         3.2          1.3         0.2  setosa
4          4.6         3.1          1.5         0.2  setosa
5          5.0         3.6          1.4         0.2  setosa
6          5.4         3.9          1.7         0.4  setosa

# head函数使用
> head(iris)
  Sepal.Length Sepal.Width Petal.Length Petal.Width Species
1          5.1         3.5          1.4         0.2  setosa
2          4.9         3.0          1.4         0.2  setosa
3          4.7         3.2          1.3         0.2  setosa
4          4.6         3.1          1.5         0.2  setosa
5          5.0         3.6          1.4         0.2  setosa
6          5.4         3.9          1.7         0.4  setosa

我们经常还会对数据进行转型处理,把数据框(data.frame)类型和其他数据类型的进行转化。我们有时会使用矩阵计算,R语言中默认供了矩阵(matrix)数据类型,可以很方便地把数据框转类型成矩阵类型,有时也需要把数据框的某一行或某一列转型为一个向量类型数据,或者把数据框变成一个list类型。通过数据的格式变换,用标准化的数据结构来满足数据分析的要求。

虽然R语言是统计语言,从性能上来说比C++/Java等语言慢不少。但对于数据分析的业务场景,用R语言来做数据处理的时候,你不用考虑程序如何架构,指针怎么定义,内存是否会泄露,只要关注你的数据和算法就行了。唯一需要注意的一点,不要直接用for循环的方式处理数据,尽量使用向量计算或矩阵计算的计算方法。当必须用循环的时候,你就需要用apply家族函数,代替for循环来做数据处理。关于apply家族函数的用法,请参考文章掌握R语言中的apply函数族

如果你的数据量比较大,1GB,10GB,甚至有100GB,对于这种规模比较大的数据集,apply的计算方式就不太能满足计算性能的要求了。你依然可以用data.table包, bigmemory包, ff包等,或者并行计算的包加速R语言在单机上的计算的性能。data.table的使用方法,请参考文章超高性能数据处理包data.table

那么再大规模的数据,超过1TB这个量级,不只是R语言,每种语言都会遇到计算性能的瓶颈。这个时候,我们需要把数据放到分布式系统中,如Hadoop或其他大数据的引擎中进行存储和计算。R语言与各种的大数据平台的通信接口都是通的,比如RHadoop,rhive, rhbase, rmongodb, rCassandra, SparkR, sparklyr等。如果你想了解hadoop的知识,请参考文章Hadoop家族系列文章RHadoop实践系列文章, R利剑NoSQL系列文章 之 Hive

2.2 时间序列类型

除了R语言的内置基础数据类型,对于金融的数据处理,一般我会把它变成标准的时间序列类型的数据,R语言中基本的时间序列的类型为 zoo 和 xts类型,当然还有一些其他包提供的数据类型。关于zoo和xts的详细介绍,请参考文章 R语言时间序列基础库zoo可扩展的时间序列xts

通过类型变换可以很方便地把的data.frame或者matrix等基础类型数据,变成xts时间序列类型的数据。时间序列类型的好处是它默认会以时间作为索引,对于量化策略来说,每条数据记录他都会有数据产生的时间,那这个时间就正好可以作为索引列的时间。


# 数据框
> df<-data.frame(A=1:6,B=rnorm(6))

# xts时间序列类型
> xdf<-xts(df,order.by=as.Date('2016-01-01')+1:6);xdf
           A           B
2016-01-02 1 -1.24013232
2016-01-03 2 -0.21014651
2016-01-04 3 -1.63251615
2016-01-05 4 -0.67279885
2016-01-06 5  0.01487863
2016-01-07 6  0.92012628

# 类型检查
> class(xdf)
[1] "xts" "zoo"

那么以时间作为数据的索引列的好处是,可以很方便地把数据以时间维度进行对齐。比如,你设计了一个股票交易策略和一个期货交易策略,由于股票是T+1交易,今天买了明天才能卖;而期货是T+0交易,今天买了马上就可以卖出。针对不同的市场规则,在设计交易策略时,可能就会选择不同的交易周期,那么这时两个策略的交易周期就会不一样,那么时间维度可能也不是对齐的。如果这两个策略是对冲的,那么我们就需要把它们以时间维度进行对齐,才能进行实现对策略模型对冲的准确计算。

把不同时间的维度的数据转化成同一个时间维度,相当于做时间的标准化。通过标准化的操作,让数据变成同一时间维度,数据之间才能够进行计算。

举个简单的例子,我们做股票交易,在实盘交易过程中,你可能最关心的是每秒最新的价格数据,每一秒都会产生一条数据,这是属于日内交易策略。另外,我们再做一个周期稍微长一点的策略,以日线为基础的,那么这里一条记录就是一天收盘价。对比日内策略,1秒钟一条数据和1天一条数据,它们不同维度的数据,是不能直接进行计算。

我们要处理这种不同周期维度数据的时候,就需要把数据转成同一个维度的。比如,我们对日线和周线的数据进行合并的时候,可以是把周线数据拆成日线数据,就是把一周分成五天。反过来,也可以把日线数据合并为周线数据,把5天的数据合并成一周。

所以这个时候就需要一个统一的数据格式进行标准化的数据定义,zoo和xts就是我们作为时间序列基础数据类型。这两个包是由第三方开发的,提供了很丰富的时间序列处理函数,我们可以直接使用这些函数来操作金融数据。很多其他的第三方金融算法分析包,也都是以这两个包作为基础开发。

3. R语言和金融模型

当我们掌握了R语言处理数据的方法,了解了如何使用R语言的基础数据类型和时间序列数据类型,下面我们就可以构建金融的策略模型。

金融建模跟其他行业的数据建模是类似的,只是由于行业不一样,金融行业有很多背景知识和金融市场规则需要我们了解。金融本身就是一个玩数据的行业,你可以通过获得交易数据,财务数据,上市公司的各种事件数据,基本面数据,宏观数据,舆情数据,互联网数据等,来构建你自己的交易策略。

我们需要把这些数据进行组合整理,结合你自己对业务的理解,使用R语言从数据中发现规律,并构建交易模型。用程序对历史数据进行回测,来验证规律的可靠性,是否会长期有效,并控制风险,最后把验证过的规律变成算法模型,这个就是金融策略建模的过程。

从金融交易分析的角度,可以从3个维度进行分析 基本面分析,技术面分析和消息面分析。

  • 基本面:指对宏观经济、行业和公司基本情况的分析,包括公司经营理念策略、公司报表等的分析。长线投资一般用基本面分析,通过基本面可以判断是否值去交易。
  • 技术面:指通过技术指标变化,判断股票走势形态,进行K线组合等,通过技术面可以判断如何进行交易。
  • 消息面:指上市公司发布的利好和利空的消息,通过消息面可以判断市场的情绪。

对于量化模型,大部分都是基于技术指标的模型,通过技术指标建模,跟踪市场的表现。在不完全了解金融业务和金融市场的情况下,通过几个技术指标来监控市场的走势,发现市场的机会也是有可能的。

量化交易和主观交易并不是对立的,量化交易是对主观交易的补充,当我们以数据作为决策基础的时候,其实可以尽量减少拍脑袋过程,创建数据模型也可以给我们心里建立良好的信心。如果交易没有使用量化的方法,那就跟我们平时做事一样,你可能想到什么就是什么。没有数据基础,那完全就是感觉,这样子交易就是很容易赔钱。

对于中国很多的散户,听到一个消息就跟着风的买卖股票,或者凭自己感觉大盘该涨了就跟进去,这些操作其实都是很不理性的。如果你通过量化的方法,即使再简单,就靠几条均线来进行判断,这样也是能给自己一个数据的基础,建立信心,而不是完全拍脑袋的事儿。

量化交易模型主要是以技术指标为主,常用的技术指标有不少,虽然简单但还是很有用的。对于很多实盘上运行的量化策略,大都会基于这些基础的指标,但并不是把每个指标单独使用。而是把多个指标通过变换组合使用,比如说MACD是均线模型,大部分的趋势策略都以MACD做为基础指标,通过变换再生成新的衍生指标。

常用的技术指标还包括KDJ、Boll、RSI、CCI等,当你直接使用这些指标的时候,可能效果并不是太好。因为市场上普遍接受了这些技术指标,已经被大量使用。单纯地用一个指标,你掌握的信息并不比别人多,所以你可能抓不到市场上赚钱的机会。

我们需要把多种技术指标或者多个维度的指标进行结合,通过组合优化的方式来降低策略的不确定风险,同时提高收益率。如果你找到了一个只有你自己知道市场规律,你的策略产生的信号完全是跟别人有区别的,你抓住了别人看不到的机会,这个才是你的赚钱机会。你领先的越多,越少人知道这个规则,那你可能赚钱的机会就越多。

建立量化模型,其实和我们平时做数据分析的思考试是一样的。要把这件事做好,我们需要把IT技术,业务知识和基础学科知识做进一步的结合,当你发现这个结合是属于你自己特有一个知识体系,你才能更好的发挥你的才能。

我们为什么要用R来做这件事情?

首先,R语言本身提供了很多数学、统计的基础包,让数学计算变得非常容易。R语言提供了常用的数据结构,向量、数据框、矩阵等,把数据变成标准化的数据,你的关注点只在数据上就可以了。另外,R语言是免费开源的,很多的第三方开发者提供了丰富的数据挖掘包,让你可以方便的使用各种算法模型,短短几行代码,就可以搞定一个复杂的事情。

R语言,在金融领域提供了很多交易框架或者计算模型,如果你了解了金融的理论知识以后,同时有一定的金融市场经验,你可以很方便的利用这些别人提供的这些技术框架,来构建自己的交易模型。CRAN上发布的金融项目,你可以去 R的官方网站 (https://cran.r-project.org/),找到Task Views 菜单里的 Finance标签。

task

通过调用第三方的程序包,自己的代码量就变的非常少。我们做一个R语言的策略,如果是很复杂的,你可能要写100-200行,但是如果你要实现同样复杂的策略,放到C++/Java去实现,这个策略就是没有1000-2000行是不可能实现的。在CRAN上面,简单数一下Finance标签下面列出的金融包就有141个,我相信没有哪种语言会比R语言对金融行业支持的更多了。

task2

虽然说R语言在性能上有些问题,但是我们会有多少了交易策略是基于一种高频的模型,对性能要求极高的呢?其实很少。就算是高频交易策略,几秒钟交易一次,R语言都可能满足要求。

海量金融数据我们怎么处理呢?

我们可以把基于海量数据的计算变成离线模型,金融行业每天都会产生大量的数据,像每日产生的交易数据,中国市场每天可能都是以GB的量来增长,跟互联网比起来不是很快,但对于你程序加载10年的数据,他要GB或TB的一个量级。

R语言本身真的很难处理这种量级的数据,但是这种量级数据对于其他语言来说同样是很难处理的。我们并不需要把这种体量的数据,都加载到内存中,进行实时数据计算。变成离线的计算模型,仅用于建模回测。把海量数据能变成离线的方式,放到hadoop或spark计算,用海量数据进行模型的训练。

我们用到的实时数据,一般就是一天或几天的数据,会不很大,每天从开盘到收盘可能也就1-2GB,对于这个大小,我们完全有能力放到内存中,进行各种各样的计算。

做量化交易难点还是在于如何发现市场机会,R语言可以很好的满足数据计算,建模,分析等的所有技术的部分。利用你的擅长,找到市场的机会,然后去实盘交易赚到钱,我们就完成了整个的交易过程。

本文并没有介绍,如何用R语言真正的去实现一个交易策略,你可以通过下面的列表找到对应的文章。

2015年我在创业,希望能推动R语言在金融量化领域的发展,但是由于种种原因项目没有持续发展。接下来,我还会以个人的方式继续努力,继续推动R在金融领域的发展。R对我们的影响和改变是非常大的,我认识R是非常好的一门语言,我会把推动R的发展,当成一项事业来做。希望也能和各位业界朋友,一起努力,把这份事业做下去。

转载请注明出处:
http://blog.fens.me/r-finance/

打赏作者