library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.4.4 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.0
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(ggthemes)
library(ggrepel)
library(xts)
## Loading required package: zoo
##
## Attaching package: 'zoo'
##
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
##
## ######################### Warning from 'xts' package ##########################
## # #
## # The dplyr lag() function breaks how base R's lag() function is supposed to #
## # work, which breaks lag(my_xts). Calls to lag(my_xts) that you type or #
## # source() into this session won't work correctly. #
## # #
## # Use stats::lag() to make sure you're not using dplyr::lag(), or you can add #
## # conflictRules('dplyr', exclude = 'lag') to your .Rprofile to stop #
## # dplyr from breaking base R's lag() function. #
## # #
## # Code in packages is not affected. It's protected by R's namespace mechanism #
## # Set `options(xts.warn_dplyr_breaks_lag = FALSE)` to suppress this warning. #
## # #
## ###############################################################################
##
## Attaching package: 'xts'
##
## The following objects are masked from 'package:dplyr':
##
## first, last
library(tsibble)
##
## Attaching package: 'tsibble'
##
## The following object is masked from 'package:zoo':
##
## index
##
## The following object is masked from 'package:lubridate':
##
## interval
##
## The following objects are masked from 'package:base':
##
## intersect, setdiff, union
theme_set(theme_minimal())
options(scipen = 6)
library(tidyverse)
library(dplyr)
chocolate <- read.csv("~/Desktop/chocolate.csv")
library(readr)
summary(chocolate)
## ref company_manufacturer company_location review_date
## Min. : 5 Length:2530 Length:2530 Min. :2006
## 1st Qu.: 802 Class :character Class :character 1st Qu.:2012
## Median :1454 Mode :character Mode :character Median :2015
## Mean :1430 Mean :2014
## 3rd Qu.:2079 3rd Qu.:2018
## Max. :2712 Max. :2021
## country_of_bean_origin specific_bean_origin_or_bar_name cocoa_percent
## Length:2530 Length:2530 Min. :0.4200
## Class :character Class :character 1st Qu.:0.7000
## Mode :character Mode :character Median :0.7000
## Mean :0.7164
## 3rd Qu.:0.7400
## Max. :1.0000
## ingredients most_memorable_characteristics rating
## Length:2530 Length:2530 Min. :1.000
## Class :character Class :character 1st Qu.:3.000
## Mode :character Mode :character Median :3.250
## Mean :3.196
## 3rd Qu.:3.500
## Max. :4.000
My chocolate csv has the review_date in the format of year, so I decided to make up a month and day and just keep the year.
chocolate$review_date <- as.Date(paste(chocolate$review_date, "-01-01", sep = ""), format = "%Y-%m-%d")
I decided to choose chocolate rating as the “response-like” variable to see how chocolate ratings have changed over the years!
average_rating <- chocolate %>%
group_by(review_date) %>%
summarize(avg_rating = mean(rating, na.rm = TRUE))
ggplot(average_rating, aes(x = review_date, y = avg_rating)) +
geom_line() +
labs(x = "Review Date", y = "Average Rating", title = "Average Chocolate Rating Over Time")
library(tsibble)
library(dplyr)
library(ggplot2)
average_ratings <- chocolate %>%
group_by(year = lubridate::year(review_date)) %>%
summarize(avg_rating = mean(rating, na.rm = TRUE))
chocolate_tsibble <- tsibble(
review_date = as.Date(paste(average_ratings$year, "-01-01", sep = ""), format = "%Y-%m-%d"),
avg_rating = average_ratings$avg_rating,
index = review_date)
ggplot(chocolate_tsibble, aes(x = review_date, y = avg_rating)) +
geom_line() +
labs(x = "Review Date", y = "Average Rating", title = "Average Chocolate Rating Over Time")
From creating this tsibble object of just the date and response variable, I ran into the issue that due to my data set only giving me the review_date in the format of year, I converted the dates to all have 01-01-year format. This then caused the issue that there are several reviews with the exact same date. Therefore, in order to resolve this I decided to take the average rating of the duplicate dates and plot that against the review_date year.
lm_model <- lm(avg_rating ~ review_date, data = chocolate_tsibble)
summary(lm_model)
##
## Call:
## lm(formula = avg_rating ~ review_date, data = chocolate_tsibble)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.13341 -0.04169 0.01386 0.02550 0.12549
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.49063611 0.16712627 14.903 0.000000000555 ***
## review_date 0.00004344 0.00001046 4.153 0.000976 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.07045 on 14 degrees of freedom
## Multiple R-squared: 0.552, Adjusted R-squared: 0.52
## F-statistic: 17.25 on 1 and 14 DF, p-value: 0.0009757
Do you need to subset the data for multiple trends?
How strong are these trends?
From this I was able to see that the coefficient for the review date is meaningful as the p-value is 0.0009757 which is really small indicating that the linear regression is something here. This suggests that there is evidence to conclude that there is a linear trend in the average rating over time. This makes much more sense as one would expect as technology becomes more advanced the chocolate becomes better!
chocolate_tsibble$avg_rating_ma <- zoo::rollmean(chocolate_tsibble$avg_rating, k = 3, fill = NA)
ggplot(chocolate_tsibble, aes(x = review_date)) +
geom_line(aes(y = avg_rating), color = "blue", alpha = 0.7) +
geom_line(aes(y = avg_rating_ma), color = "red") +
labs(x = "Review Date", y = "Average Rating", title = "Average Chocolate Rating Over Time (with Moving Average)")
## Warning: Removed 2 rows containing missing values (`geom_line()`).
I did this to try a different approach to detect seasonality in my chocolate rating over time data. I decided to calculate a rolling mean of the average rating variable in the chocolate_tsibble. The k=3 parameter specified the width or the number of observations to include in each calculation of the rolling mean, so in my case the k=3 indicates the rolling mean was calculated using the current observation and the two preceding observations. There is a difference in trend in the average rating about once every ~3 years, meaning the first 3 years there was a decrease in the average rating, then the next 3 years there was an increase, and the next 3 years there was a decrease.
library(tsibble)
library(ggplot2)
library(forecast)
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
acf_result <- Acf(chocolate_tsibble$avg_rating, lag.max = 24, plot = FALSE)
pacf_result <- Pacf(chocolate_tsibble$avg_rating, lag.max = 24, plot = FALSE)
autoplot(acf_result, main = "ACF Plot")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main`
autoplot(pacf_result, main = "PACF Plot")
## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown
## parameters: `main`
In this section I attempted to illustrate the seasonality using autocorrelation function (ACF) and partial autocorrelation function (PACF). These plots do not cross the dotted-blue lines, which indicate the bounds of the 95% confidence interval, therefore there is no evidence of any correlation with previous lags. There are no obvious, repeating patterns in either graph that indicate seasonality. This lack of seasonality is associated with the available data from the csv, the review date is listed in years which results in a smooth view of the avg_rating over time from the start. I think the lack of seasonality makes sense in this data set because we are given a general average rating of the chocolate over the years. If we had average rating on specific batches of chocolates and the dates of the cocoa beans used to produce that chocolate, then maybe we could find whether there’s a seasonality between chocolate batches (their average ratings) over time.