library(quantmod)
getSymbols(‘CPIAUCSL’, src = ‘FRED’) # CPI data
cpi <- CPIAUCSL # changing data name
Copyright By PowCoder代写 加微信 powcoder
names(cpi) <- 'index' # changing column name class(cpi) # check that class is "ts" cpi$inf <- diff(cpi$index, 12) / lag(cpi$index, 12) # inflation plot(cpi$inf) # plot inflation series - not sure if stationary cpi$inf.diff <- diff(cpi$inf) # first difference in inflation plot(cpi$inf.diff) # looks stationary data("AirPassengers") # built-in dataset in R AirPassengers # time series of year and month data of number of air passengers class(AirPassengers) plot(AirPassengers) # definitely not stationary # in the above series, three instances of nonstationary behavior # 1. overall increasing trend # 2. seasonal fluctuations # 3. size of seasonal fluctuations gets larger over time # always remember to set seed to any number before simulating random numbers set.seed(100) y <- rnorm(250) # simulating 250 standard normal random numbers y[1:10] # check the first 10 values y <- rnorm(250) # simulate again y[1:10] # check that it now looks different set.seed(100) # if you set seed back to original value y <- rnorm(250) # and simulate again y[1:10] # you'll see the original values plot(y, xlab="time", type = 'l', col = 'blue') # looks just like random noise set.seed(100) y <- 0.01 + rnorm(250, sd = 0.05) # random noise with nonzero mean and sd = 5% set.seed(100) y <- rnorm(250, mean = 0.01, sd = 0.05) # can place the mean inside the rnorm function plot(y, type = 'l', xlab = "time", col = 'blue') # also looks like random noise abline(h = 0.01, col = 'red') # horizontal line to show the mean abline(h = c(0.06, -0.04), col = 'red', lty = 2) # horizontal dashed lines to show 1 sd away from mean # Another simulation set.seed(100) y <- 1 / sqrt(3) * rt(250, df = 3) # simulating a t_3 distribution, i.e. t dist with 3 dof plot(y, xlab = "time", col = 'blue', type = 'l') # looks like random noise again but with more extreme values # Another simulation set.seed(100) y <- 1 / sqrt(3) * rt(250, df = 30) # t distribution again but with more dof plot(y, xlab = "time", col = 'blue', type = 'l') # fewer extreme values - closer to normal # Simulating a nonstationary process with linear trend set.seed(100) e <- rnorm(250) # noise values t <- 1:250 # time y <- 0.5 + 0.1 * t + e # linear trend + noise plot(t, y, col = 'blue', type = 'l', main = "Deterministic Trend + Noise") # looks... trend-y abline(a = 0.5, b = 0.1) # trendline # Simulating a Random Walk with mean 0 set.seed(100) e <- rnorm(250) # noise values y <- cumsum(e) # easy definition cbind(e, y)[1:10, ] # checking what the first 10 values look like and effect of cumsum plot(y, type = 'l', col = 'blue') # can go pretty far away from the mean; might look like a trend but it's not # another example of random walk - same as above e <- rnorm(250) y <- cumsum(e) plot(y, type = 'l', col = 'blue') # yet another example of random walk e <- rnorm(250) y <- cumsum(e) plot(y, type = 'l', col = 'blue') # correlation in an MA(1) process curve(x / (1 + x^2), -4, 4) # here, x corresponds to theta # Simulating a MA(1) Process n.obs <- 250 # number of observations mu <- 1 # mean theta <- 0.9 # parameter sigma.e <- 1 # sd of noise term set.seed(100) e <- rnorm(n.obs, sd = sigma.e) # simulate the noise terms y <- rep(0, n.obs)# Initialize y to 0 y[1] <- mu + e[1] # Starting value of y # For loop to generate the MA values: for (t in 2:n.obs) { y[t] <- mu + e[t] + theta * e[t-1] plot(y, type = 'l', col = 'blue') # hard to see any pattern # Alternate procedure instead of for loop: Vectorize! set.seed(100) # same seed e <- rnorm(n.obs, sd = sigma.e) # noise values em1 <- c(0, e[1:(n.obs - 1)]) # noise values from previous time period cbind(e, em1)[1:10,] # see what this looks like y <- mu + e + theta * em1 # generate moving average values plot(y, type = 'l', col = 'blue') # same plot! # Simulating an AR(1) process mu <- 1; sigma.e <- 1; n.obs <- 250 # same as before phi <- 0.9 # AR parameter y <- rep(0, n.obs) # initialize y set.seed(100) e <- rnorm(n.obs, sd = sigma.e) # noise values y[1] <- mu + e[1] # set y_1 # generate AR values using for loop for (t in 2:n.obs) { y[t] <- mu + phi * (y[t-1] - mu) + e[t] # Alternative to for loop set.seed(100) # same seed y <- mu + filter(e, phi, method = "recursive") # can be used for both ar and ma plot(y, type = 'l', col = 'blue') phi <- 0.99 # try different parameter y <- mu + filter(e, phi, method = "recursive") plot(y, type = 'l', col = 'blue') # should see pattern phi <- 0 # try phi = 0 y <- mu + filter(e, phi, method = "recursive") # this is just random noise plot(y, type = 'l', col = 'blue') # looks like random noise phi <- -0.99 # try negative parameter y <- mu + filter(e, phi, method = "recursive") plot(y, type = 'l', col = 'blue') # should see pattern # Function for calculating correlations ARMAacf(ar = 0.9, lag.max = 10) plot(ARMAacf(ar = 0.9, lag.max = 10), type = 'h', col = 'red') plot(ARMAacf(ar = -0.9, lag.max = 10), type = 'h', col = 'red') abline(h = 0) 程序代写 CS代考 加微信: powcoder QQ: 1823890830 Email: powcoder@163.com