第一組:黃金價格預測
報告人員:吳東翰 姚冠豪 趙致忠 段柏宇 黃丞祥 連祐瑋
2024-06-09
資料來源:https://www.kaggle.com/datasets/somyaagarwal69/gold-forecasting
目的:資料給予我們原油價格、利率（回購利率）、美元兌換盧比匯率、Sensex（BSE）、消費者物價指數和美元指數，我們去分析這些變數是否與黃金價格有相關性。
資料筆數:239筆
反應變數:Gold_Price(黃金價格)單位:千元
解釋變數:
Crude_Oil(原油價格)單位:千元
Interest_Rate(利率（回購利率）)單位:%
USD_INR(美元兌換盧比匯率)單位:盧比
Sensex(Sensex（BSE）孟買股票指數)單位:千
CPI(消費者物價指數)
USD_Index(美元指數)
讀入檔案
GoldUP <- read.csv("C:/Users/sense/大學/大二下/迴歸/GoldUP.csv")
head(GoldUP)## Date Gold_Price Crude_Oil Interest_Rate USD_INR Sensex CPI
## 1 01-10-2000 4538 1455.51 8.0 46.31830 3711.02 37.23
## 2 01-11-2000 4483 1512.47 8.0 46.78361 3997.99 37.31
## 3 01-12-2000 4541 1178.11 8.0 46.74586 3972.12 36.98
## 4 01-01-2001 4466 1208.18 8.0 46.53603 4326.72 36.90
## 5 01-02-2001 4370 1267.18 7.5 46.51459 4247.04 36.73
## 6 01-03-2001 4269 1166.45 7.0 46.60935 3604.38 36.90
## USD_Index
## 1 116.65
## 2 115.24
## 3 109.56
## 4 110.52
## 5 112.01
## 6 117.37我們沒有用到時間，因此欄位選取時沒有選擇date欄位
畫出各變數的矩陣散佈圖
data=GoldUP[,c(2:8)]
head(data)## Gold_Price Crude_Oil Interest_Rate USD_INR Sensex CPI USD_Index
## 1 4538 1455.51 8.0 46.31830 3711.02 37.23 116.65
## 2 4483 1512.47 8.0 46.78361 3997.99 37.31 115.24
## 3 4541 1178.11 8.0 46.74586 3972.12 36.98 109.56
## 4 4466 1208.18 8.0 46.53603 4326.72 36.90 110.52
## 5 4370 1267.18 7.5 46.51459 4247.04 36.73 112.01
## 6 4269 1166.45 7.0 46.60935 3604.38 36.90 117.37plot(data)
為了簡化程式碼，所以將反應變數以Y表示，解釋變數以X1~X6表示
Y=data$Gold_Price
X1=data$Crude_Oil
X2=data$Interest_Rate
X3=data$USD_INR
X4=data$Sensex
X5=data$CPI
X6=data$USD_Index所有迴歸模型比較選取法:
各種可能迴歸模型下的 \(SSE_p\),\(R^2\),\(R_a^2\),\(C_p\),\(AIC\),\(BIC\),\(PRESS_p\)
library(car)## Warning: 套件 'car' 是用 R 版本 4.3.3 來建造的## 載入需要的套件：carData## Warning: 套件 'carData' 是用 R 版本 4.3.3 來建造的source("C:/Users/sense/大學/大二下/迴歸/criteria.best.R")
library(knitr)## Warning: 套件 'knitr' 是用 R 版本 4.3.3 來建造的bestsubset(as.matrix(data[,2:7]), data$Gold_Price)## $model
## $model$R2
## [1] 64
##
## $model$R2a
## [1] 32
##
## $model$Cp
## [1] 32
##
## $model$AIC
## [1] 32
##
## $model$BIC
## [1] 32
##
##
## $stat
## p X1 X2 X3 X4 X5 X6 SSE R2 R2a Cp
## [1,] 1 0 0 0 0 0 0 32406923004 -2.220446e-16 -2.220446e-16 5803.06
## [2,] 2 1 0 0 0 0 0 18218095731 0.4378332 0.4354612 3160.521
## [3,] 2 0 1 0 0 0 0 30509111457 0.05856192 0.0545896 5451.342
## [4,] 3 1 1 0 0 0 0 17956202638 0.4459146 0.441219 3113.709
## [5,] 2 0 0 1 0 0 0 8685896612 0.731974 0.7308431 1383.893
## [6,] 3 1 0 1 0 0 0 3804333749 0.8826074 0.8816125 476.0585
## [7,] 3 0 1 1 0 0 0 8594917872 0.7347814 0.7325338 1368.936
## [8,] 4 1 1 1 0 0 0 3223131478 0.9005419 0.8992722 369.733
## [9,] 2 0 0 0 1 0 0 6303503900 0.805489 0.8046683 939.8582
## [10,] 3 1 0 0 1 0 0 5227319101 0.8386975 0.8373305 741.277
## [11,] 3 0 1 0 1 0 0 5418887618 0.8327861 0.831369 776.9819
## [12,] 4 1 1 0 1 0 0 5016717316 0.8451961 0.8432199 704.0247
## [13,] 3 0 0 1 1 0 0 5139585838 0.8414047 0.8400607 724.9251
## [14,] 4 1 0 1 1 0 0 3148500997 0.9028448 0.9016045 355.8232
## [15,] 4 0 1 1 1 0 0 4692791641 0.8551917 0.8533431 643.6508
## [16,] 5 1 1 1 1 0 0 3016700498 0.9069119 0.9053206 333.258
## [17,] 2 0 0 0 0 1 0 2599894671 0.9197735 0.919435 249.573
## [18,] 3 1 0 0 0 1 0 1916630135 0.9408574 0.9403562 124.2249
## [19,] 3 0 1 0 0 1 0 2511896607 0.9224889 0.921832 235.1717
## [20,] 4 1 1 0 0 1 0 1907734282 0.9411319 0.9403804 124.5669
## [21,] 3 0 0 1 0 1 0 2307067342 0.9288094 0.9282061 196.9953
## [22,] 4 1 0 1 0 1 0 1910578231 0.9410441 0.9402915 125.097
## [23,] 4 0 1 1 0 1 0 2181379634 0.9326879 0.9318285 175.5694
## [24,] 5 1 1 1 0 1 0 1906541214 0.9411687 0.940163 126.3445
## [25,] 3 0 0 0 1 1 0 2369414269 0.9268856 0.9262659 208.6156
## [26,] 4 1 0 0 1 1 0 1574862508 0.9514035 0.9507831 62.52567
## [27,] 4 0 1 0 1 1 0 2354630861 0.9273417 0.9264142 207.8603
## [28,] 5 1 1 0 1 1 0 1376332841 0.9575297 0.9568037 27.52336
## [29,] 4 0 0 1 1 1 0 1508651642 0.9534466 0.9528523 50.18517
## [30,] 5 1 0 1 1 1 0 1338485273 0.9586976 0.9579915 20.46926
## [31,] 5 0 1 1 1 1 0 1508594907 0.9534484 0.9526526 52.1746
## [32,] 6 1 1 1 1 1 0 1245251293 0.9615745 0.9607499 5.092146
## [33,] 2 0 0 0 0 0 1 31835493728 0.01763294 0.01348793 5698.556
## [34,] 3 1 0 0 0 0 1 14940780891 0.5389633 0.5350562 2551.689
## [35,] 3 0 1 0 0 0 1 29989588050 0.07459316 0.06675073 5356.513
## [36,] 4 1 1 0 0 0 1 13339273437 0.588382 0.5831273 2255.197
## [37,] 3 0 0 1 0 0 1 3545889234 0.8905824 0.8896551 427.8891
## [38,] 4 1 0 1 0 0 1 3015373432 0.9069528 0.905765 331.0107
## [39,] 4 0 1 1 0 0 1 3539053604 0.8907933 0.8893992 428.6151
## [40,] 5 1 1 1 0 0 1 2837761729 0.9124335 0.9109366 299.9071
## [41,] 3 0 0 0 1 0 1 6221086476 0.8080322 0.8064053 926.4971
## [42,] 4 1 0 0 1 0 1 4955375895 0.847089 0.8451369 692.5918
## [43,] 4 0 1 0 1 0 1 5345827927 0.8350406 0.8329347 765.3649
## [44,] 5 1 1 0 1 0 1 4906733201 0.84859 0.8460018 685.5256
## [45,] 4 0 0 1 1 0 1 3235201786 0.9001694 0.898895 371.9827
## [46,] 5 1 0 1 1 0 1 2780972881 0.9141858 0.9127189 289.3227
## [47,] 5 0 1 1 1 0 1 3156062716 0.9026115 0.9009467 359.2326
## [48,] 6 1 1 1 1 0 1 2734856113 0.9156089 0.9137979 282.7273
## [49,] 3 0 0 0 0 1 1 2085456710 0.9356478 0.9351024 155.6911
## [50,] 4 1 0 0 0 1 1 1882060457 0.9419241 0.9411827 119.7818
## [51,] 4 0 1 0 0 1 1 2007636811 0.9380491 0.9372583 143.1869
## [52,] 5 1 1 0 0 1 1 1882030979 0.941925 0.9409323 121.7763
## [53,] 4 0 0 1 0 1 1 2082335041 0.9357441 0.9349238 157.1093
## [54,] 5 1 0 1 0 1 1 1877658441 0.94206 0.9410695 120.9613
## [55,] 5 0 1 1 0 1 1 2007634111 0.9380492 0.9369902 145.1864
## [56,] 6 1 1 1 0 1 1 1877019437 0.9420797 0.9408368 122.8422
## [57,] 4 0 0 0 1 1 1 1546201899 0.9522879 0.9516788 57.18385
## [58,] 5 1 0 0 1 1 1 1423573225 0.9560719 0.955321 36.3281
## [59,] 5 0 1 0 1 1 1 1545555976 0.9523078 0.9514926 59.06347
## [60,] 6 1 1 0 1 1 1 1322619739 0.9591871 0.9583113 19.51221
## [61,] 5 0 0 1 1 1 1 1427726547 0.9559438 0.9551907 37.10221
## [62,] 6 1 0 1 1 1 1 1328319473 0.9590112 0.9581317 20.57454
## [63,] 6 0 1 1 1 1 1 1427209243 0.9559597 0.9550147 39.00579
## [64,] 7 1 1 1 1 1 1 1244756900 0.9615898 0.9605964 7
## AIC BIC
## [1,] 4477.317 4480.793
## [2,] 4341.663 4348.616
## [3,] 4464.894 4471.847
## [4,] 4340.202 4350.632
## [5,] 4164.632 4171.585
## [6,] 3969.323 3979.753
## [7,] 4164.116 4174.545
## [8,] 3931.7 3945.606
## [9,] 4088.01 4094.963
## [10,] 4045.267 4055.697
## [11,] 4053.87 4064.299
## [12,] 4037.439 4051.345
## [13,] 4041.222 4051.652
## [14,] 3926.101 3940.007
## [15,] 4021.486 4035.392
## [16,] 3917.881 3935.263
## [17,] 3876.343 3883.296
## [18,] 3805.472 3815.901
## [19,] 3870.114 3880.543
## [20,] 3806.36 3820.266
## [21,] 3849.784 3860.214
## [22,] 3806.716 3820.622
## [23,] 3838.396 3852.301
## [24,] 3808.21 3825.593
## [25,] 3856.157 3866.587
## [26,] 3760.532 3774.438
## [27,] 3856.662 3870.567
## [28,] 3730.328 3747.71
## [29,] 3750.266 3764.172
## [30,] 3723.663 3741.046
## [31,] 3752.257 3769.64
## [32,] 3708.407 3729.266
## [33,] 4475.065 4482.018
## [34,] 4296.264 4306.694
## [35,] 4462.789 4473.218
## [36,] 4271.166 4285.072
## [37,] 3952.509 3962.939
## [38,] 3915.776 3929.682
## [39,] 3954.048 3967.954
## [40,] 3903.266 3920.649
## [41,] 4086.864 4097.294
## [42,] 4034.499 4048.405
## [43,] 4052.625 4066.531
## [44,] 4034.141 4051.523
## [45,] 3932.594 3946.499
## [46,] 3898.435 3915.817
## [47,] 3928.674 3946.057
## [48,] 3896.439 3917.297
## [49,] 3825.648 3836.077
## [50,] 3803.121 3817.027
## [51,] 3818.559 3832.465
## [52,] 3805.118 3822.5
## [53,] 3827.29 3841.196
## [54,] 3804.562 3821.944
## [55,] 3820.558 3837.941
## [56,] 3806.48 3827.339
## [57,] 3756.142 3770.048
## [58,] 3738.393 3755.776
## [59,] 3758.042 3775.425
## [60,] 3722.814 3743.672
## [61,] 3739.09 3756.472
## [62,] 3723.841 3744.7
## [63,] 3741.003 3761.862
## [64,] 3710.312 3734.648配適最複雜模型m1，以vif檢查共線性
m1=lm(Y~X1+X2+X3+X4+X5+X6)
summary(m1)##
## Call:
## lm(formula = Y ~ X1 + X2 + X3 + X4 + X5 + X6)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3971.8 -1628.1 -197.7 1094.8 13117.3
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.545e+03 1.967e+03 -1.294 0.197060
## X1 1.277e+00 2.189e-01 5.831 1.83e-08 ***
## X2 -7.234e+02 1.833e+02 -3.946 0.000105 ***
## X3 -2.892e+02 7.590e+01 -3.809 0.000178 ***
## X4 -6.859e-01 6.319e-02 -10.856 < 2e-16 ***
## X5 6.967e+02 4.181e+01 16.665 < 2e-16 ***
## X6 -8.118e+00 2.674e+01 -0.304 0.761739
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2316 on 232 degrees of freedom
## Multiple R-squared: 0.9616, Adjusted R-squared: 0.9606
## F-statistic: 968 on 6 and 232 DF, p-value: < 2.2e-16vif(m1)## X1 X2 X3 X4 X5 X6
## 5.118266 2.104681 27.172035 21.554413 65.680179 4.059760因為在vif檢查時，有三個vif≧10的值，因此選擇將最大的值X5(vif=65.680179)刪除
以下是刪除X5後的模型配適m2模型:
m2=lm(Y~X1+X2+X3+X4+X6)
summary(m2)##
## Call:
## lm(formula = Y ~ X1 + X2 + X3 + X4 + X6)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9218.1 -2120.0 -770.4 1734.2 13140.1
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.304e+04 2.756e+03 -4.733 3.83e-06 ***
## X1 1.910e+00 3.189e-01 5.990 7.88e-09 ***
## X2 -5.364e+02 2.706e+02 -1.982 0.04864 *
## X3 7.923e+02 5.824e+01 13.603 < 2e-16 ***
## X4 1.635e-01 5.523e-02 2.961 0.00338 **
## X6 -1.790e+02 3.653e+01 -4.900 1.79e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3426 on 233 degrees of freedom
## Multiple R-squared: 0.9156, Adjusted R-squared: 0.9138
## F-statistic: 505.6 on 5 and 233 DF, p-value: < 2.2e-16vif(m2)## X1 X2 X3 X4 X6
## 4.963893 2.096794 7.312914 7.527817 3.462840檢查共線性，確定vif均＜10
運用向前選取法、向後選取法和逐步選取法，來選取最佳適配模型
none=lm(Y~1)
step(none,scope =list(upper=m2,lower=none),direction='forward' ) #向前選取法## Start: AIC=4477.32
## Y ~ 1
##
## Df Sum of Sq RSS AIC
## + X4 1 2.6103e+10 6.3035e+09 4088.0
## + X3 1 2.3721e+10 8.6859e+09 4164.6
## + X1 1 1.4189e+10 1.8218e+10 4341.7
## + X2 1 1.8978e+09 3.0509e+10 4464.9
## + X6 1 5.7143e+08 3.1835e+10 4475.1
## <none> 3.2407e+10 4477.3
##
## Step: AIC=4088.01
## Y ~ X4
##
## Df Sum of Sq RSS AIC
## + X3 1 1163918062 5139585838 4041.2
## + X1 1 1076184799 5227319101 4045.3
## + X2 1 884616281 5418887618 4053.9
## + X6 1 82417423 6221086476 4086.9
## <none> 6303503900 4088.0
##
## Step: AIC=4041.22
## Y ~ X4 + X3
##
## Df Sum of Sq RSS AIC
## + X1 1 1991084841 3148500997 3926.1
## + X6 1 1904384052 3235201786 3932.6
## + X2 1 446794197 4692791641 4021.5
## <none> 5139585838 4041.2
##
## Step: AIC=3926.1
## Y ~ X4 + X3 + X1
##
## Df Sum of Sq RSS AIC
## + X6 1 367528115 2780972881 3898.4
## + X2 1 131800499 3016700498 3917.9
## <none> 3148500997 3926.1
##
## Step: AIC=3898.44
## Y ~ X4 + X3 + X1 + X6
##
## Df Sum of Sq RSS AIC
## + X2 1 46116768 2734856113 3896.4
## <none> 2780972881 3898.4
##
## Step: AIC=3896.44
## Y ~ X4 + X3 + X1 + X6 + X2##
## Call:
## lm(formula = Y ~ X4 + X3 + X1 + X6 + X2)
##
## Coefficients:
## (Intercept) X4 X3 X1 X6 X2
## -1.304e+04 1.635e-01 7.923e+02 1.910e+00 -1.790e+02 -5.364e+02step(m2,scope =list(upper=m2,lower=none),direction='backward' ) #向後選取法## Start: AIC=3896.44
## Y ~ X1 + X2 + X3 + X4 + X6
##
## Df Sum of Sq RSS AIC
## <none> 2734856113 3896.4
## - X2 1 46116768 2780972881 3898.4
## - X4 1 102905616 2837761729 3903.3
## - X6 1 281844384 3016700498 3917.9
## - X1 1 421206603 3156062716 3928.7
## - X3 1 2171877088 4906733201 4034.1##
## Call:
## lm(formula = Y ~ X1 + X2 + X3 + X4 + X6)
##
## Coefficients:
## (Intercept) X1 X2 X3 X4 X6
## -1.304e+04 1.910e+00 -5.364e+02 7.923e+02 1.635e-01 -1.790e+02step(none,scope =list(upper=m2,lower=none),direction='both' ) #逐步選取法## Start: AIC=4477.32
## Y ~ 1
##
## Df Sum of Sq RSS AIC
## + X4 1 2.6103e+10 6.3035e+09 4088.0
## + X3 1 2.3721e+10 8.6859e+09 4164.6
## + X1 1 1.4189e+10 1.8218e+10 4341.7
## + X2 1 1.8978e+09 3.0509e+10 4464.9
## + X6 1 5.7143e+08 3.1835e+10 4475.1
## <none> 3.2407e+10 4477.3
##
## Step: AIC=4088.01
## Y ~ X4
##
## Df Sum of Sq RSS AIC
## + X3 1 1.1639e+09 5.1396e+09 4041.2
## + X1 1 1.0762e+09 5.2273e+09 4045.3
## + X2 1 8.8462e+08 5.4189e+09 4053.9
## + X6 1 8.2417e+07 6.2211e+09 4086.9
## <none> 6.3035e+09 4088.0
## - X4 1 2.6103e+10 3.2407e+10 4477.3
##
## Step: AIC=4041.22
## Y ~ X4 + X3
##
## Df Sum of Sq RSS AIC
## + X1 1 1991084841 3148500997 3926.1
## + X6 1 1904384052 3235201786 3932.6
## + X2 1 446794197 4692791641 4021.5
## <none> 5139585838 4041.2
## - X3 1 1163918062 6303503900 4088.0
## - X4 1 3546310774 8685896612 4164.6
##
## Step: AIC=3926.1
## Y ~ X4 + X3 + X1
##
## Df Sum of Sq RSS AIC
## + X6 1 367528115 2780972881 3898.4
## + X2 1 131800499 3016700498 3917.9
## <none> 3148500997 3926.1
## - X4 1 655832753 3804333749 3969.3
## - X1 1 1991084841 5139585838 4041.2
## - X3 1 2078818104 5227319101 4045.3
##
## Step: AIC=3898.44
## Y ~ X4 + X3 + X1 + X6
##
## Df Sum of Sq RSS AIC
## + X2 1 46116768 2734856113 3896.4
## <none> 2780972881 3898.4
## - X4 1 234400551 3015373432 3915.8
## - X6 1 367528115 3148500997 3926.1
## - X1 1 454228905 3235201786 3932.6
## - X3 1 2174403014 4955375895 4034.5
##
## Step: AIC=3896.44
## Y ~ X4 + X3 + X1 + X6 + X2
##
## Df Sum of Sq RSS AIC
## <none> 2734856113 3896.4
## - X2 1 46116768 2780972881 3898.4
## - X4 1 102905616 2837761729 3903.3
## - X6 1 281844384 3016700498 3917.9
## - X1 1 421206603 3156062716 3928.7
## - X3 1 2171877088 4906733201 4034.1##
## Call:
## lm(formula = Y ~ X4 + X3 + X1 + X6 + X2)
##
## Coefficients:
## (Intercept) X4 X3 X1 X6 X2
## -1.304e+04 1.635e-01 7.923e+02 1.910e+00 -1.790e+02 -5.364e+02★發現最佳適配模型就是m2
以下是殘差分析: 畫關係圖
par(mfrow=c(2,2))
plot(m2)
1.檢查是否具均齊性
\(H_0\):變異數具均齊性
\(H_1\):變異數不具均齊性
ncvTest(m2)## Non-constant Variance Score Test
## Variance formula: ~ fitted.values
## Chisquare = 56.09481, Df = 1, p = 6.9059e-14library(lmtest)## Warning: 套件 'lmtest' 是用 R 版本 4.3.3 來建造的## 載入需要的套件：zoo## Warning: 套件 'zoo' 是用 R 版本 4.3.3 來建造的##
## 載入套件：'zoo'## 下列物件被遮斷自 'package:base':
##
## as.Date, as.Date.numericbptest(m2)##
## studentized Breusch-Pagan test
##
## data: m2
## BP = 59.275, df = 5, p-value = 1.716e-11►ncvTest 檢定結果之p值<0.05，推翻虛無假設，表示變異數不具均齊性
►bptest 檢定結果之p值<0.05，推翻虛無假設，表示變異數不具均齊性
2.檢定殘差是否具常態性
\(H_0\):殘差具常態性
\(H_1\):殘差不具常態性
e=m2$residuals #殘差
qqnorm(e)
qqline(e)
shapiro.test(e)##
## Shapiro-Wilk normality test
##
## data: e
## W = 0.95131, p-value = 3.441e-07►因為shapiro.test結果之p值<0.05，所以推翻虛無假設，表示殘差不具常態性
3.檢定殘差是否具一階自我相關
\(H_0\):殘差無一階自我相關
\(H_1\):殘差有一階自我相關
durbinWatsonTest(m2) #雙尾檢定## lag Autocorrelation D-W Statistic p-value
## 1 0.915927 0.1024024 0
## Alternative hypothesis: rho != 0dwtest(m2) #單尾檢定##
## Durbin-Watson test
##
## data: m2
## DW = 0.1024, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0►由於兩種檢定結果之p值<0.05，推翻虛無假設，表示殘差具一階自我相關
★因為不具均齊性、常態性，所以做Box-Cox轉換
做Box-Cox轉換
trans=powerTransform(m2)
summary(trans)## bcPower Transformation to Normality
## Est Power Rounded Pwr Wald Lwr Bnd Wald Upr Bnd
## Y1 0.3051 0.33 0.1401 0.4701
##
## Likelihood ratio test that transformation parameter is equal to 0
## (log transformation)
## LRT df pval
## LR test, lambda = (0) 12.40575 1 0.00042801
##
## Likelihood ratio test that no transformation is needed
## LRT df pval
## LR test, lambda = (1) 72.32213 1 < 2.22e-16\(\lambda\)=0.3050862，以\(\lambda\)=0.3取代來做newy
newy=bcPower(Y,0.3)
newy## [1] 38.34766 38.19547 38.35593 38.14816 37.87862 37.59053 37.58477 38.07836
## [9] 37.96329 37.90689 38.10072 38.60210 38.75629 38.52861 38.45207 38.76975
## [17] 39.32110 39.37064 39.68304 40.15067 40.36659 40.05558 39.90695 40.18310
## [25] 40.11069 40.19306 40.68709 41.41990 41.46420 40.66038 40.06561 40.96880
## [33] 40.84411 40.48956 40.73070 41.34037 41.28639 41.64729 42.20210 42.39172
## [41] 42.03099 41.95849 41.79894 41.38252 41.67495 42.12348 42.27370 42.37173
## [49] 42.79167 43.20074 42.97356 42.32500 42.23346 42.57512 42.32946 42.05810
## [57] 42.29156 42.12123 42.54651 43.16661 43.87749 44.47440 45.29424 45.93621
## [65] 46.14592 46.18467 47.82709 49.44725 47.77065 48.78644 48.63873 47.84929
## [73] 47.32602 48.09028 48.07846 47.97012 48.76352 48.47512 48.39369 47.64355
## [81] 47.34698 47.40627 47.57108 48.37704 49.00132 50.02901 49.98407 51.45631
## [89] 52.30978 53.33238 52.20000 52.66519 52.95397 53.85956 52.27183 52.77148
## [97] 53.41165 52.66519 53.72088 54.43875 56.09003 56.60516 55.69552 55.87350
## [105] 55.89536 55.99349 56.27247 57.17837 57.36130 58.67493 58.78594 58.28933
## [113] 58.08828 58.18733 58.27494 59.70926 60.45546 59.96617 60.15985 60.77113
## [121] 61.18978 61.84400 62.19710 61.90607 62.01918 62.51076 63.10647 63.71271
## [129] 63.78896 64.16093 66.93473 68.21787 67.57151 69.04689 68.66442 68.30443
## [137] 68.70750 68.50314 69.04689 69.24996 69.98738 69.72739 70.25255 71.32624
## [145] 70.88691 71.23772 70.72663 70.50028 70.09342 69.76158 68.41714 67.65196
## [153] 67.97561 67.77819 70.27887 70.46615 70.63787 70.65954 70.05002 69.81206
## [161] 70.27522 69.99697 69.62383 69.17964 68.17098 68.62206 68.67751 67.80590
## [169] 67.82806 67.12717 67.57948 68.15220 67.89443 67.17004 67.58347 67.90863
## [177] 67.55395 66.65730 66.81310 67.23302 67.49883 66.74678 66.38309 67.03232
## [185] 68.80955 69.22200 69.24014 69.85432 69.93274 70.80508 71.03998 70.97426
## [193] 70.17272 69.97041 68.42567 69.18569 69.57603 69.19477 69.38486 68.97687
## [201] 69.24392 68.74667 69.29675 70.04633 69.75489 69.67602 69.19856 70.01319
## [209] 70.40652 70.42762 70.83669 70.92990 70.63570 70.19178 69.85802 70.51334
## [217] 71.19586 70.85392 71.02642 71.71773 72.39996 71.58192 71.31634 71.37217
## [225] 72.33423 73.28056 75.11557 75.47339 75.65183 75.59473 75.57609 76.73202
## [233] 77.45192 78.08733 79.89930 80.35675 80.87952 82.22286 83.75446以newy對單一變數來看是否顯著
summary(lm(newy~X1))##
## Call:
## lm(formula = newy ~ X1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -18.741 -6.311 -2.769 6.220 32.819
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.729e+01 1.387e+00 26.88 <2e-16 ***
## X1 6.110e-03 3.715e-04 16.45 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.892 on 237 degrees of freedom
## Multiple R-squared: 0.533, Adjusted R-squared: 0.531
## F-statistic: 270.5 on 1 and 237 DF, p-value: < 2.2e-16summary(lm(newy~X2))##
## Call:
## lm(formula = newy ~ X2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -23.748 -12.360 2.776 8.169 33.103
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 37.9069 4.6550 8.143 2.19e-14 ***
## X2 2.9987 0.6826 4.393 1.68e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 12.51 on 237 degrees of freedom
## Multiple R-squared: 0.07531, Adjusted R-squared: 0.0714
## F-statistic: 19.3 on 1 and 237 DF, p-value: 1.685e-05summary(lm(newy~X3))##
## Call:
## lm(formula = newy ~ X3)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.6946 -5.2550 -0.0741 4.7285 17.4734
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.77601 2.68312 1.407 0.161
## X3 1.00864 0.04898 20.593 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 7.791 on 237 degrees of freedom
## Multiple R-squared: 0.6415, Adjusted R-squared: 0.64
## F-statistic: 424.1 on 1 and 237 DF, p-value: < 2.2e-16summary(lm(newy~X4))##
## Call:
## lm(formula = newy ~ X4)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.812 -4.106 -1.839 2.852 13.502
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.876e+01 7.033e-01 55.12 <2e-16 ***
## X4 1.061e-03 3.310e-05 32.05 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5.634 on 237 degrees of freedom
## Multiple R-squared: 0.8126, Adjusted R-squared: 0.8118
## F-statistic: 1027 on 1 and 237 DF, p-value: < 2.2e-16summary(lm(newy~X6))##
## Call:
## lm(formula = newy ~ X6)
##
## Residuals:
## Min 1Q Median 3Q Max
## -18.0089 -12.1883 0.1215 11.4094 26.2691
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 83.8642 6.5624 12.779 < 2e-16 ***
## X6 -0.2863 0.0722 -3.965 9.72e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 12.6 on 237 degrees of freedom
## Multiple R-squared: 0.06221, Adjusted R-squared: 0.05825
## F-statistic: 15.72 on 1 and 237 DF, p-value: 9.723e-05►發現各變數都具顯著性
以newy做反應變數，X1、X2、X3、X4、X6做解釋變數配適m3
m3=lm(newy~X1+X2+X3+X4+X6)
summary(m3)##
## Call:
## lm(formula = newy ~ X1 + X2 + X3 + X4 + X6)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.6177 -2.9213 -0.4184 2.3459 9.5242
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.395e+01 2.797e+00 12.141 < 2e-16 ***
## X1 1.734e-03 3.236e-04 5.358 2.02e-07 ***
## X2 1.996e-01 2.746e-01 0.727 0.468
## X3 7.019e-01 5.910e-02 11.875 < 2e-16 ***
## X4 3.419e-04 5.605e-05 6.100 4.38e-09 ***
## X6 -3.007e-01 3.707e-02 -8.111 2.86e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.477 on 233 degrees of freedom
## Multiple R-squared: 0.9298, Adjusted R-squared: 0.9283
## F-statistic: 617.4 on 5 and 233 DF, p-value: < 2.2e-16畫關係圖
par(mfrow=c(2,2))
plot(m3)
1.檢查是否具均齊性
\(H_0\):變異數具均齊性
\(H_1\):變異數不具均齊性
library(car)
ncvTest(m3)## Non-constant Variance Score Test
## Variance formula: ~ fitted.values
## Chisquare = 14.51749, Df = 1, p = 0.00013886library(lmtest)
bptest(m3)##
## studentized Breusch-Pagan test
##
## data: m3
## BP = 60.367, df = 5, p-value = 1.021e-11►ncvTest 檢定結果之p值<0.05，推翻虛無假設，表示變異數不具均齊性
►bptest 檢定結果之p值<0.05，推翻虛無假設，表示變異數不具均齊性
2.檢定殘差是否具常態性
\(H_0\):殘差具常態性
\(H_1\):殘差不具常態性
e1=m3$residuals
qqnorm(e1)
qqline(e1)
shapiro.test(e1)##
## Shapiro-Wilk normality test
##
## data: e1
## W = 0.97522, p-value = 0.0003411►因為shapiro.test結果之p值<0.05，所以推翻虛無假設，表示殘差不具常態性
3.檢定殘差是否具一階自我相關
\(H_0\):殘差無一階自我相關
\(H_1\):殘差有一階自我相關
durbinWatsonTest(m3)## lag Autocorrelation D-W Statistic p-value
## 1 0.9421082 0.1042253 0
## Alternative hypothesis: rho != 0dwtest(m3)##
## Durbin-Watson test
##
## data: m3
## DW = 0.10423, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0►由於兩種檢定結果之p值均<0.05，推翻虛無假設，表示殘差具一階自我相關
★因為m3模型X2變數不顯著，殘差不具均齊和常態性，因此刪除X2變數配適m4模型
m4 = lm(newy ~ X1+X3+X4+X6)
summary(m4)##
## Call:
## lm(formula = newy ~ X1 + X3 + X4 + X6)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.8349 -2.9274 -0.4099 2.4485 9.2902
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.383e+01 2.789e+00 12.132 < 2e-16 ***
## X1 1.888e-03 2.436e-04 7.753 2.75e-13 ***
## X3 7.143e-01 5.652e-02 12.637 < 2e-16 ***
## X4 3.223e-04 4.906e-05 6.568 3.26e-10 ***
## X6 -2.938e-01 3.579e-02 -8.207 1.52e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.473 on 234 degrees of freedom
## Multiple R-squared: 0.9297, Adjusted R-squared: 0.9285
## F-statistic: 773.2 on 4 and 234 DF, p-value: < 2.2e-16畫關係圖
par(mfrow=c(2,2))
plot(m4)
1.檢查是否具均齊性
\(H_0\):變異數具均齊性
\(H_1\):變異數不具均齊性
ncvTest(m4)## Non-constant Variance Score Test
## Variance formula: ~ fitted.values
## Chisquare = 15.069, Df = 1, p = 0.00010365bptest(m4)##
## studentized Breusch-Pagan test
##
## data: m4
## BP = 40.68, df = 4, p-value = 3.131e-08►ncvTest 檢定結果之p值<0.05，推翻虛無假設，表示變異數不具均齊性
►bptest 檢定結果之p值<0.05，推翻虛無假設，表示變異數不具均齊性
2.檢定殘差是否具常態性
\(H_0\):殘差具常態性
\(H_1\):殘差不具常態性
qqnorm(m4$resid)
qqline(m4$resid)
shapiro.test(m4$resid)##
## Shapiro-Wilk normality test
##
## data: m4$resid
## W = 0.97605, p-value = 0.0004536►因為shapiro.test結果之p值<0.05，所以推翻虛無假設，表示殘差不具常態性
3.檢定殘差是否具一階自我相關
\(H_0\):殘差無一階自我相關
\(H_1\):殘差有一階自我相關
durbinWatsonTest(m4) ## lag Autocorrelation D-W Statistic p-value
## 1 0.9422079 0.1051106 0
## Alternative hypothesis: rho != 0dwtest(m4)##
## Durbin-Watson test
##
## data: m4
## DW = 0.10511, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0►由於兩種檢定結果之p值均<0.05，推翻虛無假設，表示殘差具一階自我相關
★藉由m3跟m4的比較，我們發現兩者差異不大(均不具均齊性和常態性，並且具一階自我相關)，所以選擇m3模型當最終模型，因為多一個解釋變數X2
所以我們使用m3做離群值、影響點檢測
p=5 #解釋變數個數
n=nrow(data)
ei = m3$resid
hii = hatvalues(m3)
ri = rstandard(m3)
ti = rstudent(m3)
Di = cooks.distance(m3)
DFFITSi = dffits(m3)
DFBETASi = dfbetas(m3)
case1 = which(hii>2*(p+1)/n); names(case1)=NULL; case1## [1] 1 8 9 155 156 172 234 235 236 237 238 239case2 = which(abs(ri)>2); names(case2)=NULL; case2## [1] 131 132 134 215 216 217 235 236case3 = which(abs(ti)>qt(1-0.05/n, n-p-1)); names(case3)=NULL; case3## integer(0)case4 = which(Di>qf(0.5, p+1, n-(p+1))); names(case4)=NULL; case4## integer(0)case5 = which(abs(DFFITSi)>1); names(case5)=NULL; case5## integer(0)case6 = which(abs(DFBETASi)>1, arr.ind=T); names(case6)=NULL; case6## row col我們看到case1和case2有重複的點(第235和236筆資料)，因此我們先刪除第235筆資料，重新配適m3_1
m3_1=lm(newy[-235]~X1[-235]+X2[-235]+X3[-235]+X4[-235]+X6[-235])
summary(m3_1)##
## Call:
## lm(formula = newy[-235] ~ X1[-235] + X2[-235] + X3[-235] + X4[-235] +
## X6[-235])
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.5221 -2.9012 -0.4152 2.4091 9.4754
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.356e+01 2.781e+00 12.065 < 2e-16 ***
## X1[-235] 1.802e-03 3.227e-04 5.584 6.56e-08 ***
## X2[-235] 2.410e-01 2.732e-01 0.882 0.379
## X3[-235] 6.806e-01 5.949e-02 11.440 < 2e-16 ***
## X4[-235] 3.495e-04 5.573e-05 6.272 1.73e-09 ***
## X6[-235] -2.911e-01 3.706e-02 -7.856 1.47e-13 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.45 on 232 degrees of freedom
## Multiple R-squared: 0.9304, Adjusted R-squared: 0.9289
## F-statistic: 619.8 on 5 and 232 DF, p-value: < 2.2e-16ncvTest(m3_1)## Non-constant Variance Score Test
## Variance formula: ~ fitted.values
## Chisquare = 15.05292, Df = 1, p = 0.00010454bptest(m3_1)##
## studentized Breusch-Pagan test
##
## data: m3_1
## BP = 63.739, df = 5, p-value = 2.047e-12qqnorm(m3_1$residuals)
qqline(m3_1$residuals)
shapiro.test(m3_1$resid)##
## Shapiro-Wilk normality test
##
## data: m3_1$resid
## W = 0.9767, p-value = 0.0005882durbinWatsonTest(m3_1) ## lag Autocorrelation D-W Statistic p-value
## 1 0.93983 0.1070578 0
## Alternative hypothesis: rho != 0dwtest(m3_1)##
## Durbin-Watson test
##
## data: m3_1
## DW = 0.10706, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0★因為殘差分析出來的結果一樣不具均齊性、常態性且具一階自我相關，所以我們還是選用m3當最終模型
\(\hat{GoldPrice}\)=33.9544 +0.0017\(CrudeOil\)+0.1996\(InterestRate\)+0.7019\(USDINR\)+3^{-4}\(Sensex\)+-0.3007\(USDIndex\)