# ==============================================================================
# VALENCIA DANA: NARRATIVES, EMOTIONS, AND ENGAGEMENT ON X
# Analysis Script
# ==============================================================================


# 1. LOAD PACKAGES----------------------------------------------------------

library(readr)
library(dplyr)
library(AER)
library(MASS)


# 2. DATA IMPORT-------------------------------------------

DANA_sentim_etiquetas_autor_X <- read_csv(
  "DANA NOV25/DANA_sentim_etiquetas_autor_X.csv"
)


# 3. VARIABLE CONSTRUCTION-------------------------------------------------

# Engagement Rate

DANA_sentim_etiquetas_autor_X <- DANA_sentim_etiquetas_autor_X %>%
  mutate(
    eng_rate = (
      retweetcount +
      replycount +
      likecount +
      quotecount
    ) / viewcount
  )


# Identify the Dominant Narrative

DANA_sentim_etiquetas_autor_X <- DANA_sentim_etiquetas_autor_X %>%
  mutate(
    narrativa = case_when(
      pred_bulos_desinfo == 1 ~ "Desinformación / Bulos",
      pred_clima_causa == 1 ~ "Causa climática",
      pred_blame_pol == 1 ~ "Culpabilización política",
      pred_ayuda_solidaridad == 1 ~ "Ayuda / Solidaridad",
      pred_info_util == 1 ~ "Información útil",
      pred_otros_no_clasif == 1 ~ "Otros / No clasificado",
      TRUE ~ "Sin narrativa identificada"
    )
  )

DANA_sentim_etiquetas_autor_X %>%
  count(narrativa) %>%
  mutate(
    porcentaje = 100 * n / sum(n)
  )
# Define the Event Date

fecha_evento <- as.Date("2024-10-29")


# Create the Variable Measuring Days Since the Event

DANA_sentim_etiquetas_autor_X <- DANA_sentim_etiquetas_autor_X %>%
  mutate(
    fecha_tuit = as.Date(createdat),
    dias_evento = as.numeric(
      difftime(
        fecha_tuit,
        fecha_evento,
        units = "days"
      )
    )
  )


# Express the Number of Followers in Hundreds of Thousands

DANA_sentim_etiquetas_autor_X <- DANA_sentim_etiquetas_autor_X %>%
  mutate(
    followers_k = author_followers / 100000
  )

mean(DANA_sentim_etiquetas_autor_X$author_followers, na.rm = TRUE)

# Avoid Scientific Notation in the Output
options(scipen = 999)


# Convert the Narrative Variable to a Categorical Variable
DANA_sentim_etiquetas_autor_X$narrativa <- factor(
  DANA_sentim_etiquetas_autor_X$narrativa
)


# ==============================================================================
# 4. NARRATIVE MODELS
# ==============================================================================


# 4.1. OLS MODEL ---------------------------------------------------------------

modelo_ols <- lm(
  eng_rate ~
    dias_evento +
    narrativa,
  data = DANA_sentim_etiquetas_autor_X,
  contrasts = list(
    narrativa = contr.treatment(
      levels(DANA_sentim_etiquetas_autor_X$narrativa),
      base = which(
        levels(DANA_sentim_etiquetas_autor_X$narrativa) ==
          "Información útil"
      )
    )
  )
)


# Results

summary(modelo_ols)


# Goodness-of-Fit Measures

summary(modelo_ols)$r.squared
summary(modelo_ols)$adj.r.squared

AIC(modelo_ols)
BIC(modelo_ols)


# RMSE

sqrt(mean(modelo_ols$residuals^2))


# Set "Información Útil" as the Reference Category

DANA_sentim_etiquetas_autor_X$narrativa <- relevel(
  DANA_sentim_etiquetas_autor_X$narrativa,
  ref = "Información útil"
)


# 4.2. Tobit model -------------------------------------------------------------

modelo_tobit <- tobit(
  eng_rate ~
    dias_evento +
    narrativa,
  left = 0,
  data = DANA_sentim_etiquetas_autor_X
)


# Results

summary(modelo_tobit)


# Log-likelihood and McFadden pseudo-R² 

LL <- as.numeric(logLik(modelo_tobit))
LL0 <- modelo_tobit$null.logLik

R2_mcfadden_tobit <- 1 - (LL / LL0)

R2_mcfadden_tobit
logLik(modelo_tobit)


# Goodness-of-Fit Measures

AIC(modelo_tobit)
BIC(modelo_tobit)


# 4.3. Descriptive Analysis Prior to the Ordered Logit Model ----------------------------

summary(DANA_sentim_etiquetas_autor_X$eng_rate)

quantile(
  DANA_sentim_etiquetas_autor_X$eng_rate,
  probs = seq(0, 1, 0.1),
  na.rm = TRUE
)

hist(
  DANA_sentim_etiquetas_autor_X$eng_rate,
  breaks = 50
)

quantile(
  DANA_sentim_etiquetas_autor_X$eng_rate,
  probs = c(0, 0.25, 0.50, 0.75, 0.90, 0.95, 0.99, 1),
  na.rm = TRUE
)


# # Create Ordinal Engagement Categories

DANA_sentim_etiquetas_autor_X$eng_cat <- cut(
  DANA_sentim_etiquetas_autor_X$eng_rate,
  breaks = c(
    -Inf,
    0.006378003,
    0.047718176,
    Inf
  ),
  labels = c(
    "bajo",
    "medio",
    "alto"
  ),
  include.lowest = TRUE
)


# Convert to an Ordered Factor

DANA_sentim_etiquetas_autor_X$eng_cat <- ordered(
  DANA_sentim_etiquetas_autor_X$eng_cat,
  levels = c(
    "bajo",
    "medio",
    "alto"
  )
)


# 4.4. Ordered logit model -----------------------------------------------------

modelo_ologit <- polr(
  eng_cat ~
    dias_evento +
    narrativa,
  data = DANA_sentim_etiquetas_autor_X,
  Hess = TRUE
)


# Results

summary(modelo_ologit)


# Log-likelihood 

LL <- as.numeric(logLik(modelo_ologit))


# Null model

modelo_null <- update(
  modelo_ologit,
  . ~ 1
)

LL0 <- as.numeric(logLik(modelo_null))


# McFadden Pseudo-R² 

R2_mcfadden_ologit <- 1 - (LL / LL0)

R2_mcfadden_ologit


# Goodness-of-Fit Measures

AIC(modelo_ologit)
BIC(modelo_ologit)
logLik(modelo_ologit)


# ==============================================================================
# 5.EMOTION MODELS
# ==============================================================================


# 5.1. OLS Model ---------------------------------------------------------------

modelo01 <- lm(
  eng_rate ~
    dias_evento +
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    nlp_sadness_score,
  data = DANA_sentim_etiquetas_autor_X
)


# Results

summary(modelo01)
AIC(modelo01)


# 5.2. Tobit model -------------------------------------------------------------

modelo_tobit_01 <- tobit(
  eng_rate ~
    dias_evento +
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score,
  left = 0,
  data = DANA_sentim_etiquetas_autor_X
)


# Results

summary(modelo_tobit_01)
AIC(modelo_tobit_01)


# 5.3. Ordered logit model-----------------------------------------------------

modelo_ologit_01 <- polr(
  eng_cat ~
    dias_evento +
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score,
  data = DANA_sentim_etiquetas_autor_X,
  Hess = TRUE
)


# Results

summary(modelo_ologit_01)
logLik(modelo_ologit_01)


# P-values

s <- summary(modelo_ologit_01)

pvals <- 2 * (
  1 - pnorm(
    abs(s$coefficients[, "t value"])
  )
)


# Results

resultado <- cbind(
  s$coefficients,
  "p-value" = pvals
)

resultado




# ==============================================================================
# 6. MODELS WITH NARRATIVES AND EMOTIONS
# ==============================================================================


# 6.1. OLS Model ---------------------------------------------------------------

modelo2 <- lm(
  eng_rate ~
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    nlp_sadness_score +
    dias_evento +
    narrativa,
  data = DANA_sentim_etiquetas_autor_X
)


# Results

summary(modelo2)
AIC(modelo2)


# 7.2. Tobit model -------------------------------------------------------------

modelo_tobit_02 <- tobit(
  eng_rate ~
    nlp_anger_score +
    dias_evento +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    narrativa,
  left = 0,
  data = DANA_sentim_etiquetas_autor_X
)


# Results

summary(modelo_tobit_02)
AIC(modelo_tobit_02)


# 7.3. Ordered logit model-----------------------------------------------------

modelo_ologit_02 <- polr(
  eng_cat ~
    dias_evento +
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    narrativa,
  data = DANA_sentim_etiquetas_autor_X,
  Hess = TRUE
)


# Results

summary(modelo_ologit_02)


# P-values

s <- summary(modelo_ologit_02)

pvals <- 2 * (
  1 - pnorm(
    abs(s$coefficients[, "t value"])
  )
)


# Results

resultado <- cbind(
  s$coefficients,
  "p-value" = pvals
)

resultado


# Goodness-of-Fit Measures

AIC(modelo_ologit_02)
logLik(modelo_ologit_02)


# ==============================================================================
# 7. MODELS WITH NARRATIVE–TIME INTERACTIONS
# ==============================================================================


# 7.1. OLS model ---------------------------------------------------------------

modelo_full <- lm(
  eng_rate ~
    dias_evento +
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    narrativa +
    narrativa * dias_evento,
  data = DANA_sentim_etiquetas_autor_X
)


# Results

summary(modelo_full)
AIC(modelo_full)


# 7.2. Tobit model -------------------------------------------------------------

modelo_tobit_03 <- tobit(
  eng_rate ~
    dias_evento +
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    narrativa +
    narrativa * dias_evento,
  left = 0,
  data = DANA_sentim_etiquetas_autor_X
)


# Results

summary(modelo_tobit_03)
AIC(modelo_tobit_03)


# 7.3. Ordered logit model-----------------------------------------------------

modelo_ologit_03 <- polr(
  eng_cat ~
    dias_evento +
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    narrativa +
    narrativa * dias_evento,
  data = DANA_sentim_etiquetas_autor_X,
  Hess = TRUE
)


# Results

summary(modelo_ologit_03)


# P-values

s <- summary(modelo_ologit_03)

pvals <- 2 * (
  1 - pnorm(
    abs(s$coefficients[, "t value"])
  )
)


# Results

resultado <- cbind(
  s$coefficients,
  "p-value" = pvals
)

resultado


# # Goodness-of-Fit Measures

AIC(modelo_ologit_03)
logLik(modelo_ologit_03)


# ==============================================================================
# 8. NEGATIVE BINOMIAL MODELS FOR RETWEET COUNTS
# ==============================================================================


# 8.1. Narrative model ---------------------------------------------------------

modelo_rt_nb1 <- glm.nb(
  retweetcount ~
    narrativa +
    dias_evento +
    offset(log(viewcount + 1)),
  data = DANA_sentim_etiquetas_autor_X
)

# Results

summary(modelo_rt_nb1)

# Goodness-of-Fit Measures

AIC(modelo_rt_nb1)
BIC(modelo_rt_nb1)
logLik(modelo_rt_nb1)



# 8.2. Emotion model -----------------------------------------------------------

modelo_rt_nb2 <- glm.nb(
  retweetcount ~
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    dias_evento +
    offset(log(viewcount + 1)),
  data = DANA_sentim_etiquetas_autor_X
)

# Results

summary(modelo_rt_nb2)

# Goodness-of-Fit Measures

AIC(modelo_rt_nb2)
BIC(modelo_rt_nb2)
logLik(modelo_rt_nb2)



# 8.3. Narratives and emotions model -------------------------------------------

modelo_rt_nb3 <- glm.nb(
  retweetcount ~
    narrativa +
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    dias_evento +
    offset(log(viewcount + 1)),
  data = DANA_sentim_etiquetas_autor_X
)

# Results

summary(modelo_rt_nb3)

# Goodness-of-Fit Measures

AIC(modelo_rt_nb3)
BIC(modelo_rt_nb3)
logLik(modelo_rt_nb3)



# 8.4. Narrative–time interaction model ----------------------------------------

modelo_rt_nb4 <- glm.nb(
  retweetcount ~
    narrativa +
    nlp_anger_score +
    nlp_joy_score +
    nlp_optimism_score +
    nlp_hate_score +
    dias_evento +
    narrativa * dias_evento +
    offset(log(viewcount + 1)),
  data = DANA_sentim_etiquetas_autor_X
)

# Results

summary(modelo_rt_nb4)

# Goodness-of-Fit Measures

AIC(modelo_rt_nb4)
BIC(modelo_rt_nb4)
logLik(modelo_rt_nb4)