3
$\begingroup$

I have a dataset of 2 variables that should be heavily correlated.

There are some underlying reasons why this set has an R^2 of only 0.620 when modeled in a simple Linear Regression; the independent variable data has 'bled' into itself by the way it has been produced.

What I'd like is to see if I can use autocorrelation to help find some errors I can propagate into the independent variable to resolve some of these issues but being quite new to these methods I'm unsure of how to implement them.

I've tried a moving 2-day average on the independent variable that I feel is causing the issue. This did increase the R^2 which adjusted to about 0.787 after averaged but I'd like to do much better than that.

I've been using this Python script below to help me understand autocorrelation better.

Here is the dataset I've been working on including the 2-day averaging: Test Set for West Data

What I understand so far is that autocorrelation can find a lag term to fit the least error between a set using self-similarity. What I need to figure out is how to find these errors for the length of my independent variable.

Below is my first attempt at using this script to lag my series but because it's a ML algorithm the prediction errors calculated are only the length of the test set. I'm sure it's due to my lack of understanding but how can I use these errors and propagate them into my original length 67 independent variable set to help fit my model better?

enter image description here

Above is the graph of the independent variable's predicted values (red) using autocorrelation on the test set (blue).

$\endgroup$
2
$\begingroup$

Since the data has a seasonality of 5, we used seasonality of 5.

Here is a plot of the Y and X showing a strong positive relationship and some visually obvious outliers.enter image description here

The ACF/PACF looks as follows:

enter image description here

enter image description here

Here is the forecast.

enter image description here

Here is the model. The model includes the causal plus 4 seasonal dummies and 11 outliers and a change in the seasonality for periods 2 and 3 that went down severely beginning at period 37 and 58.

enter image description here

Here are the residual ACF/PACF

enter image description here

enter image description here

Here is the dataset cleansed of the outliers

[![enter image description here][9]][9]
| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ Is there a way you could explain each piece of what you're doing? As much as I like seeing an answer the reason I posted my question was to understand better how lag terms work and how to apply error terms to better fit seasonal models. $\endgroup$ – HelloToEarth Nov 20 '18 at 23:52
  • $\begingroup$ Sorry, but I had the Y and X backwards....so I reposted....our system recognizes seasonality of 5 as daily data and considers seasonal dummies automatically. See the other answer without the seasonal dummies included. $\endgroup$ – Tom Reilly Nov 21 '18 at 15:33
1
$\begingroup$

If I may add some more info on top of this discussion, I would like to show what follows.

test_ds <- read.csv("test.csv", header = TRUE, stringsAsFactors = FALSE)
test_ds$Date <- as.Date(test_ds$Date, format=c("%Y-%m-%d"))
head(test_ds)

        Date West.Counts West.Demands  X X2.Day.West.Counts.Avg
1 2005-01-01         452        67797 NA                    NaN
2 2005-02-01         333        56151 NA                  392.5
3 2005-03-01         262        40801 NA                  297.5
4 2005-11-01         222        18694 NA                  242.0
5 2005-12-01         523        48155 NA                  372.5
6 2006-01-01         339        52230 NA                  431.0

The cross-correlation plot shows what follows ((I do not think prewhitening is needed):

ccf(x=test_ds$West.Counts, y=test_ds$West.Demands)

enter image description here

Then cross-correlation at lag = -1 may suggest an additional regressor. Also we can try to add some dummy seasonal variable.

wc_ts <- ts(test_ds$West.Counts, frequency=5)
ll <- nrow(test_ds)-1
west_counts_lag <- c(NA, test_ds$West.Counts[1:ll])
dummy_seas <- factor(cycle(wc_ts))

lm_res <- lm(West.Demands ~ West.Counts + west_counts_lag + dummy_seas, data= test_ds)
summary(lm_res)



Call:
lm(formula = West.Demands ~ West.Counts + west_counts_lag + dummy_seas, 
    data = test_ds)

Residuals:
   Min     1Q Median     3Q    Max 
-40817  -3166    535   3322  13837 

Coefficients:
                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)      13437.78    7098.42   1.893 0.063175 .  
West.Counts         62.57      12.31   5.084 3.88e-06 ***
west_counts_lag     41.72      12.24   3.407 0.001176 ** 
dummy_seas2      -3471.68    3570.92  -0.972 0.334851    
dummy_seas3      -4192.77    4525.39  -0.926 0.357899    
dummy_seas4     -16984.64    4601.65  -3.691 0.000484 ***
dummy_seas5      -2614.11    3886.59  -0.673 0.503785    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 7349 on 60 degrees of freedom
  (1 observation deleted due to missingness)
Multiple R-squared:  0.8566,    Adjusted R-squared:  0.8422 
F-statistic: 59.71 on 6 and 60 DF,  p-value: < 2.2e-16

Diagnostic plots:

enter image description here

enter image description here

enter image description here

enter image description here

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ You do nice work. If your error plot vs time shows any outliers, then the model you are estimating is not robust. $\endgroup$ – Tom Reilly Nov 21 '18 at 15:20
0
$\begingroup$

Since the data is frequency of 5, our default is daily data so I reran suppressing the automatic search for day of the week dummies.

We now have a lag of the X variable. 9 outliers were identified and a level shift higher beginning at period 38 and forward.

enter image description here

| cite | improve this answer | |
$\endgroup$
0
$\begingroup$

Forcing intercept equals zero (i.e. "0 +" in the regression formula) helps in reaching higher R^2. Moreover, West.Demands acf plot may suggest a regressor at lag = -2.

enter image description here

test_ds <- read.csv("test.csv", header = TRUE, stringsAsFactors = FALSE)
test_ds$Date <- as.Date(test_ds$Date, format=c("%Y-%m-%d"))

ll <- nrow(test_ds)-1
west_counts_lag <- c(NA, test_ds$West.Counts[1:ll])
ll2 <- ll - 1
west_demands_lag2 <- c(NA, NA, test_ds$West.Demands[1:ll2])

lm_res <- lm(West.Demands ~ 0 + West.Counts + west_counts_lag + west_demands_lag2, data= test_ds)
summary(lm_res)

Call:
lm(formula = West.Demands ~ 0 + West.Counts + west_counts_lag + 
    west_demands_lag2, data = test_ds)

Residuals:
   Min     1Q Median     3Q    Max 
-42281  -2333   1324   4915  13886 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
West.Counts       75.39004    6.17358  12.212  < 2e-16 ***
west_counts_lag   59.09928    7.51237   7.867  6.1e-11 ***
west_demands_lag2 -0.07648    0.03980  -1.922   0.0591 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 8295 on 63 degrees of freedom
  (2 observations deleted due to missingness)
Multiple R-squared:  0.9747,    Adjusted R-squared:  0.9735 
F-statistic: 809.8 on 3 and 63 DF,  p-value: < 2.2e-16

Diagnostic plots:

enter image description here enter image description here

enter image description here enter image description here

| cite | improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.