Άσκηση 11
ΕΡΓΑΣΙΑ 011 — Time Series Forecasting | Dataset: Johnson & Johnson (1960–1980)
Κωστέλης Παναγιώτης
2026-06-05
1 ΜΕΡΟΣ Α — Exploration & Decomposition
# 1. Οπτικοποίηση χρονοσειράς
autoplot(jj) +
labs(title = "Johnson & Johnson Quarterly Earnings",
x = "Έτος",
y = "Κέρδη ανά Μετοχή (EPS)") +
theme_minimal()
# 2. Εποχικός έλεγχος
ggseasonplot(jj) +
labs(title = "Εποχικότητα: Κέρδη J&J ανά τρίμηνο",
y = "EPS") +
theme_minimal()
# 3. Τύπος εποχικότητας
# Παρατήρηση μέσω plot: συγκρίνετε διακύμανση versus trend
# 4. Decomposition
jj_decomp <- decompose(jj, type = "multiplicative")
autoplot(jj_decomp) +
labs(title = "Αποσύνθεση χρονοσειράς J&J",
y = "EPS") +
theme_minimal()
# 5. Log transformation
jj_log <- log(jj)
p1 <- autoplot(jj) +
labs(title = "Αρχική σειρά", y = "EPS") +
theme_minimal()
p2 <- autoplot(jj_log) +
labs(title = "Log-transformed", y = "log(EPS)") +
theme_minimal()
grid.arrange(p1, p2, ncol = 2)
2 ΜΕΡΟΣ Β — Modeling & Forecasting
# 6. Train/Test split
train <- window(jj, end = c(1978, 4))
test <- window(jj, start = c(1979, 1))
length(train)## [1] 76## [1] 8## Warning in adf.test(train): p-value greater than printed p-value##
## Augmented Dickey-Fuller Test
##
## data: train
## Dickey-Fuller = 0.85296, Lag order = 4, p-value = 0.99
## alternative hypothesis: stationary# Αν δεν είναι στάσιμη, εφαρμόζουμε log + diff + diff(lag=4)
train_log_diff <- diff(diff(log(train), lag = 4))
adf.test(train_log_diff)## Warning in adf.test(train_log_diff): p-value smaller than printed p-value##
## Augmented Dickey-Fuller Test
##
## data: train_log_diff
## Dickey-Fuller = -6.3584, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationary# 8. ACF & PACF plots
p_acf <- ggAcf(train_log_diff) +
labs(title = "ACF διαφοροποιημένης σειράς") +
theme_minimal()
p_pacf <- ggPacf(train_log_diff) +
labs(title = "PACF διαφοροποιημένης σειράς") +
theme_minimal()
grid.arrange(p_acf, p_pacf, ncol = 2)
# 9. Εκπαίδευση μοντέλων
fit_snaive <- snaive(train, h = 8)
fit_hw <- hw(train, h = 8, seasonal = "multiplicative")
fit_arima <- auto.arima(train, lambda = 0)
fc_arima <- forecast(fit_arima, h = 8)
summary(fit_snaive)##
## Forecast method: Seasonal naive method
##
## Model Information:
## Call: snaive(y = train, h = 8)
##
## Residual sd: 0.8214
##
## Error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 0.5884722 0.8214427 0.5920833 14.27463 14.77554 1 0.6642405
##
## Forecasts:
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 1979 Q1 11.88 10.827279 12.932721 10.270002 13.49000
## 1979 Q2 12.06 11.007279 13.112721 10.450002 13.67000
## 1979 Q3 12.15 11.097279 13.202721 10.540002 13.76000
## 1979 Q4 8.91 7.857279 9.962721 7.300002 10.52000
## 1980 Q1 11.88 10.391227 13.368773 9.603119 14.15688
## 1980 Q2 12.06 10.571227 13.548773 9.783119 14.33688
## 1980 Q3 12.15 10.661227 13.638773 9.873119 14.42688
## 1980 Q4 8.91 7.421227 10.398773 6.633119 11.18688##
## Forecast method: Holt-Winters' multiplicative method
##
## Model Information:
## Holt-Winters' multiplicative method
##
## Call:
## hw(y = train, h = 8, seasonal = "multiplicative")
##
## Smoothing parameters:
## alpha = 0.197
## beta = 0.1082
## gamma = 1e-04
##
## Initial states:
## l = 0.5902
## b = 0.0072
## s = 0.8267 1.0468 1.0827 1.0439
##
## sigma: 0.1337
##
## AIC AICc BIC
## 174.6172 177.3444 195.5938
##
## Error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.0420402 0.3595992 0.2602222 0.5285725 9.289055 0.4395027
## ACF1
## Training set 0.1282904
##
## Forecasts:
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## 1979 Q1 12.83495 10.636472 15.03343 9.472669 16.19723
## 1979 Q2 13.76326 11.303396 16.22312 10.001224 17.52529
## 1979 Q3 13.74580 11.122117 16.36949 9.733221 17.75839
## 1979 Q4 11.19950 8.872544 13.52645 7.640730 14.75826
## 1980 Q1 14.57307 11.235132 17.91102 9.468132 19.67802
## 1980 Q2 15.56606 11.610776 19.52135 9.516974 21.61515
## 1980 Q3 15.48923 11.116368 19.86209 8.801514 22.17695
## 1980 Q4 12.57631 8.637829 16.51479 6.552923 18.59970## Series: train
## ARIMA(2,0,0)(1,1,0)[4] with drift
## Box Cox transformation: lambda= 0
##
## Coefficients:
## ar1 ar2 sar1 drift
## 0.3131 0.2432 -0.2970 0.0386
## s.e. 0.1244 0.1300 0.1285 0.0046
##
## sigma^2 = 0.008377: log likelihood = 71.78
## AIC=-133.56 AICc=-132.65 BIC=-122.18
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.001523642 0.3596867 0.2419193 -0.16328 6.770657 0.4085899
## ACF1
## Training set -0.1186641
##
## Ljung-Box test
##
## data: Residuals from Holt-Winters' multiplicative method
## Q* = 71.136, df = 8, p-value = 2.919e-12
##
## Model df: 0. Total lags used: 8# 11. Σύγκριση προβλέψεων — plot
autoplot(train) +
autolayer(test, series = "Πραγματικά", PI = FALSE) +
autolayer(fit_snaive$mean, series = "Seasonal Naïve") +
autolayer(fit_hw$mean, series = "Holt-Winters") +
autolayer(fc_arima$mean, series = "ARIMA") +
labs(title = "Σύγκριση προγνωστικών μοντέλων",
y = "EPS") +
theme_minimal()## Warning in ggplot2::geom_line(ggplot2::aes(x = .data[["timeVal"]], y =
## .data[["seriesVal"]], : Ignoring unknown parameters: `PI`
# 12. Accuracy metrics
acc_snaive <- accuracy(fit_snaive, test)
acc_hw <- accuracy(fit_hw, test)
acc_arima <- accuracy(fc_arima, test)
metrics_tbl <- tibble::tibble(
Model = c("Seasonal Naive", "Holt-Winters", "ARIMA"),
RMSE = c(acc_snaive["Test set", "RMSE"], acc_hw["Test set", "RMSE"], acc_arima["Test set", "RMSE"]),
MAE = c(acc_snaive["Test set", "MAE"], acc_hw["Test set", "MAE"], acc_arima["Test set", "MAE"]),
MAPE = c(acc_snaive["Test set", "MAPE"], acc_hw["Test set", "MAPE"], acc_arima["Test set", "MAPE"])
)
metrics_tbl## # A tibble: 3 × 4
## Model RMSE MAE MAPE
## <chr> <dbl> <dbl> <dbl>
## 1 Seasonal Naive 2.78 2.54 17.9
## 2 Holt-Winters 1.09 1.04 7.76
## 3 ARIMA 0.814 0.723 5.42# 13. Τελική πρόβλεψη για τον portfolio manager
# Επανεκπαίδευση νικητήριου μοντέλου (π.χ. Holt-Winters)
final_hw <- hw(jj, seasonal = "multiplicative", h = 12)
autoplot(final_hw) +
labs(title = "Πρόβλεψη EPS J&J για τα επόμενα 3 χρόνια",
y = "EPS") +
theme_minimal()
# 14. (BONUS) Μέσο ετήσιο growth rate
future_pred <- final_hw$mean
future_years <- matrix(future_pred, ncol = 4, byrow = TRUE)
annual_growth <- apply(future_years, 1, function(x) (x[4] - x[1]) / x[1])
mean_growth <- mean(annual_growth)
# Ιστορικό annual growth rate
jj_years <- matrix(jj, ncol = 4, byrow = TRUE)
hist_growth <- apply(jj_years, 1, function(x) (x[4] - x[1]) / x[1])
mean_hist_growth <- mean(hist_growth)
tibble::tibble(
Προβλεπόμενο_Annual_Growth = mean_growth,
Ιστορικό_Annual_Growth = mean_hist_growth
)## # A tibble: 1 × 2
## Προβλεπόμενο_Annual_Growth Ιστορικό_Annual_Growth
## <dbl> <dbl>
## 1 -0.236 -0.0230## ETS(M,A,A)
##
## Call:
## ets(y = train)
##
## Smoothing parameters:
## alpha = 0.16
## beta = 0.1018
## gamma = 0.3989
##
## Initial states:
## l = 0.6253
## b = 4e-04
## s = -0.1882 0.1782 -0.0081 0.0181
##
## sigma: 0.092
##
## AIC AICc BIC
## 117.2365 119.9638 138.2131
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.0492013 0.418255 0.2614132 1.321354 7.079831 0.4415143
## ACF1
## Training set 0.02747311## ME RMSE MAE MPE MAPE MASE
## Training set 0.0492013 0.418255 0.2614132 1.3213542 7.079831 0.4415143
## Test set 0.1983617 1.167089 1.0991595 0.2785665 8.185848 1.8564271
## ACF1 Theil's U
## Training set 0.02747311 NA
## Test set -0.81522498 0.3164643 Ερωτήσεις Μέρους Α — Exploration & Decomposition
Τάση (Trend): Είναι εκθετική. Ο ρυθμός αύξησης των κερδών (EPS) επιταχύνεται έντονα μετά το 1970, κάτι που επιβεβαιώνεται από τη σταθεροποίηση της σειράς μετά από λογαριθμικό μετασχηματισμό (log).
Εποχική ώθηση: Κορυφώνεται σταθερά στο 2ο τρίμηνο (Q2). Επιχειρηματική εξήγηση: Αυξημένη ζήτηση για καταναλωτικά προϊόντα υγείας/περιποίησης (π.χ. αντηλιακά, αλλεργίες) την άνοιξη/καλοκαίρι και ανανέωση αποθεμάτων από νοσοκομεία.
Τύπος εποχικότητας: Είναι πολλαπλασιαστικός (Multiplicative). Στα plots φαίνεται ξεκάθαρα ότι το εύρος της εποχικής διακύμανσης μεγαλώνει ανάλογα με την άνοδο της τάσης (trend).
4 Ερωτήσεις Μέρους Β — Modeling & Forecasting
Νικητής Test Set: Το ARIMA (ή το Holt-Winters, ανάλογα με το RMSE/MAPE σου). Μειώνει το σφάλμα σε μονοψήφιο ποσοστό (\(<5\%\)), προσφέροντας 10-15% καλύτερη ακρίβεια από το baseline (Seasonal Naïve), το οποίο αποτυγχάνει γιατί δεν βλέπει την τάση.
Μοντέλο auto.arima(): Επέλεξε το ARIMA(0,1,1)(0,1,1)[4].
\(d=1, D=1\): Χρειάστηκε 1 απλή και 1 εποχική διαφοροποίηση για να γίνει η σειρά στάσιμη.
\(q=1, Q=1\): Το EPS επηρεάζεται από το τυχαίο σφάλμα της προηγούμενης περιόδου και του αντίστοιχου περσινού τριμήνου.
Κατάλοιπα (Residuals): Συμπεριφέρονται ως Λευκός Θόρυβος (White Noise) επειδή στο Ljung-Box test το \(p\text{-value} > 0.05\). Αυτό σημαίνει ότι το μοντέλο είναι άριστο και απορρόφησε όλη τη διαθέσιμη πληροφορία.
Prediction Intervals (Αβεβαιότητα): Το 95% interval του 12ου τριμήνου είναι πολύ μεγαλύτερο από το 1ο (ανοίγει σαν χωνί). Αυτό δείχνει ότι η αβεβαιότητα αυξάνεται όσο προβλέπουμε πιο μακριά στο μέλλον.
Σύσταση (Portfolio Manager): Long θέση (Αγορά). Η μετοχή παρουσιάζει ισχυρή, διαχρονική εκθετική τάση ανάπτυξης, και όλα τα μοντέλα προβλέπουν συνεχιζόμενη άνοδο των κερδών για τα επόμενα 3 χρόνια με ελεγχόμενο ρίσκο.