Replace a specific column in the data structure.

Prepare updated data in "w". v_XCPG_1992_2016 only has data from 1992. Use "w["1992::"]"  to limit data just from 1992. Data were updated according to GDPnow 2016/5/28.


> tail(v_XCPG_1992_2016$GDP)
2015-11-01 16477.960
2015-12-01 16485.354
2016-01-01 16505.700
2016-02-01 16543.000
2016-03-01 16576.000
2016-04-01 16607.000
> tail(w)
2015-11-01 16482.07
2015-12-01 16493.57
2016-01-01 16505.07
2016-02-01 16540.61
2016-03-01 16580.61
2016-04-01 16620.61
> v_XCPG_1992_2016$GDP <- w["1992::"]
> tail(v_XCPG_1992_2016$GDP)
                GDP
2015-11-01 16482.07
2015-12-01 16493.57
2016-01-01 16505.07
2016-02-01 16540.61
2016-03-01 16580.61
2016-04-01 16620.61

S&P500 forecast next 10 months

With newly updated GDP data and model at  May 25th , S&P500 forecast next 10 months are as below.

Please note that the current price is already more than 100 pts below theoretical value. 2323.292 is expected  at the end of this year.


> 0.603*(predict(VAR(v_XCPG_1992_2016,lag=6))$fcst$GDP[,1])-7791
 [1] 2238.212 2254.216 2268.441 2282.354 2293.338 2303.527 2313.299 2323.292 2332.697 2341.712

GDP forecast by newly updated data.

2016 1st Quarter GDP was revised upward. The growth rate was updated from 0.5% to 0.8%.
The new GDP forecast are as below. The previous one is available at May 23rd entry. #10 months is up 6.1 trillion #8 is also up 6.7. It will also give S&P 500 from another 30 to 40 pts at the end of this year.



> predict(VAR(v_XCPG_1992_2016,lag=6))$fcst$GDP
          fcst    lower    upper        CI
 [1,] 16632.19 16602.29 16662.10  29.90217
 [2,] 16658.73 16593.67 16723.80  65.06516
 [3,] 16682.32 16577.21 16787.44 105.11452
 [4,] 16705.40 16566.14 16844.65 139.25646
 [5,] 16723.61 16552.99 16894.24 170.62578
 [6,] 16740.51 16541.47 16939.54 199.03576
 [7,] 16756.71 16529.33 16984.10 227.38612
 [8,] 16773.29 16518.24 17028.34 255.05082
 [9,] 16788.88 16506.56 17071.21 282.32406
[10,] 16803.83 16495.07 17112.60 308.76632

check normal distribution (shapilo.test)

Outlook often gives viewers a false impression. Shapiro.test function will test whether data are according to a normal distribution. When p-Value is more than 0.05, the distribution can be called "normal".

Samples below shows when data are extracted just from Monday to Friday without considering national holiday, it's not a normal distribution. However, other groups 1) workday (business day only), 2) all-day and 3) weekend (Saturday and Sunday) distribution are normal. The degree of each distribution is different, though.


> shapiro.test(i[weekday])

Shapiro-Wilk normality test
data:  i[weekday]
W = 0.97156, p-value = 0.02788

> shapiro.test(i[workday])

Shapiro-Wilk normality test
data:  i[workday]
W = 0.98811, p-value = 0.5692

> shapiro.test(i)

Shapiro-Wilk normality test

data:  i

W = 0.98489, p-value = 0.1182

> shapiro.test(i[weekend])

Shapiro-Wilk normality test

data:  i[weekend]

W = 0.95236, p-value = 0.07851

do not confuse with "plot.mlm is not implemented yet"

Please note that lm output class differs according to a given formula. In most cases, it just belongs to "lm". However, it seems that when the target parameter has multiple culumns, result output also belongs to "mlm"

> class(lm(to.monthly(SP500["2009::2016-02-01"]) ~ GDP_XTS["2009::2016-02-01"]))
[1] "mlm" "lm" 
> plot(lm(to.monthly(SP500["2009::2016-02-01"]) ~ GDP_XTS["2009::2016-02-01"]))
 エラー: 'plot.mlm' is not implemented yet
> plot(lm(to.monthly(SP500["2009::2016-02-01"])[,6] ~ GDP_XTS["2009::2016-02-01"]))
> class(lm(to.monthly(SP500["2009::2016-02-01"])[,6] ~ GDP_XTS["2009::2016-02-01"]))
[1] "lm"
> plot(lm(to.monthly(SP500["2009::2016-02-01"])[,6] ~ GDP_XTS["2009::2016-02-01"]))
>

lm(to.quarterly(SP5["2010-04-01::2016-03-01"])[,6] ~ GDPC96["2010-04-01::2016-01-01"] )

A strong correlation exists on those periods, however, the degree of correlation is lower in different terms. It seems that the relationship between GDP and  stock price is more likely to become solid in an upward trend.

> summary(lm(to.quarterly(SP5["2010-04-01::2016-03-01"])[,4] ~ GDPC96["2010-04-01::2016-01-01"] ) )

Call:
lm(formula = to.quarterly(SP5["2010-04-01::2016-03-01"])[, 4] ~  GDPC96["2010-04-01::2016-01-01"])

Residuals:
    Min      1Q  Median      3Q     Max 
-196.63  -85.43   12.53   71.48  171.59 

Coefficients:
                                                      Estimate       Std. Error     t value Pr(>|t|)    
(Intercept)                                       -7.791e+03  6.168e+02  -12.63 1.47e-11 ***
GDPC96["2010-04-01::2016-01-01"]  6.036e-01  3.954e-02   15.27 3.45e-13 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 105.7 on 22 degrees of freedom
Multiple R-squared:  0.9137, Adjusted R-squared:  0.9098 
F-statistic:   233 on 1 and 22 DF,  p-value: 3.453e-13

This model indicates that S&P 500 will up 248 points when GDP annual growth rate is 2.5%.

> last(as.numeric(GDPC96)) * 0.025 * 0.603

[1] 248.6282

Fear for low birth rate and hyper inflation and Gereman balanced budget policy

Fear for low birth rate and hyper inflation and Gereman balanced budget policy. Previously thought that German long-going and unchallenged balanced policy came from hyper inflation after WW1, but now need necessity to see the relationship with low birthrate at those periods. It makes a dramatic contrast with Japan's failure to maintain the policy after "Oil Shock and stagnation" while they didn't have experienced low birthrate in inter-war  period.

calculate X months after the specific date.

"seq" function included in mondate library is enhanced to calculate the distance of months between two dates. This is not possible in the standard "seq" function. Please don't forget to load mondate library beforehand.


> seq(mondate("2009/12/31"), by=1, length.out=8)[8]
mondate: timeunits="months"
[1] 2010/07/31

predict by VAR analysis

In some cases, as VAR predict function's output  all implied data on the console, it may take extra time to find the data of focus. Specifying attribute as below suppresses unnecessary outputs.

>predict(VAR(v_XCPG_1992_2016),lag=6)$fcst$<tag of focus data>

Below will help to understand output's data structure.

>predict(VAR(v_XCPG_1992_2016),lag=6)[1]
>predict(VAR(v_XCPG_1992_2016),lag=6)[2]
>predict(VAR(v_XCPG_1992_2016),lag=6)[3]

Please see the sample outputs.

> predict(VAR(v_XCPG_1992_2016),lag=6)$fcst$GDP
          fcst         lower         upper        CI
 [1,] 16577.50 16526.29 16628.71  51.21148
 [2,] 16596.18 16523.71 16668.65  72.46879
 [3,] 16614.76 16525.95 16703.58  88.81411
 [4,] 16633.27 16530.64 16735.89 102.62429
 [5,] 16651.68 16536.86 16766.50 114.81930
 [6,] 16670.02 16544.15 16795.89 125.87053
 [7,] 16688.28 16552.22 16824.34 136.05765
 [8,] 16706.46 16560.90 16852.03 145.56265
 [9,] 16724.58 16570.07 16879.09 154.51178
[10,] 16742.63 16579.63 16905.63 162.99681
> 

update fundamental datas 統計 主要統計データ FRED YAHOO 指標

update data by getsymbols.

PAYEMS is All Employees: Total Nonfarm Payrolls
RRSFS is Real Retail and Food Services Sales
UNDCONTSA is New Privately-Owned Housing Units Under Construction: Total
PI is Personal Income
TWEXB is Trade Weighted U.S. Dollar Index: Broad

each operation will give data in XTS format.

getSymbols("GDPC96",src="FRED")
getSymbols("PAYEMS",src="FRED")
getSymbols("UNDCONTSA",src="FRED");
getSymbols("RRSFS",src="FRED")
getSymbols("PI",src="FRED")
getSymbols("TWEXB",src="FRED")

# PCE

getSymbols("PCEPI",src="FRED",auto.assign=TRUE)

# CPI

getSymbols("CPIAUCSL",src="FRED",auto.assign=TRUE)

# Core CPI

getSymbols("CPILFESL",src="FRED",auto.assign=TRUE)

# core PCE

getSymbols("PCEPILFE",src="FRED",auto.assign=TRUE)

how to convert from CSV to XTS

when time series data is given in CSV format as below.

~$ head SP5.csv 

Date,Open,High,Low,Close,Volume,Adj Close

2016-05-02,2067.169922,2084.870117,2025.910034,2052.320068,4068314600,2052.320068

2016-04-01,2056.620117,2111.050049,2033.800049,2065.300049,4087129000,2065.300049

2016-03-01,1937.089966,2072.209961,1937.089966,2059.73999,4379759000,2059.73999

the combination of read.csv, read.zoo and as.xts is able to transform and conver t data into xts objects.

> SP5 <- as.xts(read.zoo(read.csv("~/SP5.csv")))

> tail(SP5)

              Open    High     Low   Close     Volume Adj.Close

2015-12-01 2082.93 2104.27 1993.26 2043.94 3922935900   2043.94

2016-01-04 2038.20 2038.20 1812.29 1940.24 5153017800   1940.24

2016-02-01 1936.94 1962.96 1810.10 1932.23 4881887000   1932.23

2016-03-01 1937.09 2072.21 1937.09 2059.74 4379759000   2059.74

2016-04-01 2056.62 2111.05 2033.80 2065.30 4087129000   2065.30

2016-05-02 2067.17 2084.87 2025.91 2052.32 4068314600   2052.32

> class(SP5)

[1] "xts" "zoo"

plot() function and frequency conversion function like to.monthly() are available as usual.

> plot(to.yearly(SP5))