Three motivating questions

These examples illustrate the kinds of research questions sentiment analysis can address in environmental science:

1. Tone over time. Has US newspaper coverage of wildfire become more negative over the past two decades? If so, does the shift correlate with increases in burned acreage, or does media framing lead or lag the physical trend?

2. Entity-level framing. Within a single article about the Inflation Reduction Act, is the tone toward the EPA positive while the tone toward the fossil-fuel industry is negative? A single “article-level” sentiment score would miss this entirely.

3. Cross-outlet comparison. Do articles from a regional paper in a coal-dependent community frame energy transition differently than articles from a national outlet? Comparing sentiment distributions across sources reveals structural differences in coverage.

What sentiment analysis is not

Sentiment analysis does not tell you what a reader thinks, feels, or believes. It tells you something about the language in the text itself, as measured by a particular dictionary or model. Conflating text-level polarity with audience reception is a common analytical error.

It also does not handle framing well. Two articles can have identical sentiment scores while framing the same issue in completely different ways (economic growth vs. ecological collapse). Frame analysis requires different tools, though sentiment can serve as one input.

When to use which method

The decision depends on your sample size, research question, and computational resources:

This tutorial follows the pipeline above. We start with the simplest methods (lexicons) and add complexity only when needed.

Simulated environmental news corpus

For this tutorial, we create a small simulated corpus of environmental news headlines and lead paragraphs. In a real project, you would import text from files, a database, or an API.

# ── Simulated environmental news articles ──────────────────────────────
articles <- tibble(
  doc_id = paste0("art_", sprintf("%03d", 1:20)),
  source = rep(c("National Tribune", "Valley Gazette",
                 "Coastal Herald", "Mountain Post"), each = 5),
  date   = seq(as.Date("2024-01-15"), by = "month", length.out = 20),
  topic  = rep(c("wildfire", "water", "energy", "biodiversity", "climate"),
               times = 4),
  text   = c(
    # National Tribune
    "Devastating wildfire season scorches 2.3 million acres across the West, forcing thousands of evacuations and overwhelming firefighting resources.",
    "New EPA water quality standards will impose stricter limits on agricultural runoff, drawing praise from conservation groups and criticism from farm lobbies.",
    "Solar installations surged 42% in 2023, marking a record year for clean energy deployment and creating thousands of new jobs in underserved communities.",
    "Alarming decline in migratory bird populations raises concerns about habitat loss and pesticide exposure across major flyways.",
    "Global temperatures shattered records for the twelfth consecutive month, intensifying calls for accelerated emissions reductions at the upcoming climate summit.",
    # Valley Gazette
    "Local ranchers struggle to balance wildfire risk with grazing traditions as prescribed burn programs expand into private lands.",
    "Community groups celebrate the restoration of riparian habitat along Bear Creek after years of advocacy and volunteer effort.",
    "The proposed wind farm faces opposition from residents concerned about noise, visual impact, and declining property values in the rural township.",
    "Native plant restoration project transforms abandoned lot into thriving pollinator garden, attracting dozens of butterfly species.",
    "Farmers report shorter growing seasons and unpredictable frost dates, attributing the changes to shifting climate patterns.",
    # Coastal Herald
    "After years of devastating fire seasons, a coalition of agencies unveils a bold new strategy for community wildfire preparedness and resilient landscapes.",
    "Toxic algal bloom shuts down water intake for coastal city, exposing vulnerabilities in aging infrastructure and nutrient management.",
    "Offshore wind project receives final permits after environmental review finds minimal impact on marine mammals and seabirds.",
    "Kelp forest recovery surprises scientists as urchin removal program shows rapid regrowth of canopy and return of marine biodiversity.",
    "Sea level projections released today paint a stark picture for low-lying neighborhoods, with managed retreat now part of the planning conversation.",
    # Mountain Post
    "Controlled burns reduce fuel loads but spark air quality warnings in mountain communities already burdened by summer smoke.",
    "Drought emergency declared as reservoir levels drop below 30% capacity, triggering mandatory conservation measures for the first time in a decade.",
    "Utility announces early retirement of its last coal plant, pledging to replace capacity with battery storage and community solar.",
    "Endangered species listing for mountain caribou prompts heated debate between conservation advocates and logging industry representatives.",
    "Glaciers in the range have lost 40% of their volume since 1980, according to a university survey released this week."
  )
)

# Build a quanteda corpus with document-level variables
corp <- corpus(articles, text_field = "text",
               docid_field = "doc_id")

# Quick look at the corpus
summary(corp, n = 5)

#> Corpus consisting of 20 documents, showing 5 documents:
#> 
#>     Text Types Tokens Sentences           source       date        topic
#>  art_001    20     20         1 National Tribune 2024-01-15     wildfire
#>  art_002    23     24         1 National Tribune 2024-02-15        water
#>  art_003    25     26         1 National Tribune 2024-03-15       energy
#>  art_004    18     18         1 National Tribune 2024-04-15 biodiversity
#>  art_005    20     22         1 National Tribune 2024-05-15      climate

Tokenization and document-feature matrix

Most quanteda operations require tokenizing the corpus (splitting text into words) and converting it to a document-feature matrix (dfm), where rows are documents and columns are unique word types.

# Tokenize: lowercase, remove punctuation and stopwords
toks <- tokens(corp,
               remove_punct  = TRUE,
               remove_numbers = TRUE) %>%
  tokens_tolower() %>%
  tokens_remove(stopwords("en"))

# Document-feature matrix
dfm_articles <- dfm(toks)
dfm_articles

#> Document-feature matrix of: 20 documents, 253 features (94.39% sparse) and 3 docvars.
#>          features
#> docs      devastating wildfire season scorches million acres across west
#>   art_001           1        1      1        1       1     1      1    1
#>   art_002           0        0      0        0       0     0      0    0
#>   art_003           0        0      0        0       0     0      0    0
#>   art_004           0        0      0        0       0     0      1    0
#>   art_005           0        0      0        0       0     0      0    0
#>   art_006           0        1      0        0       0     0      0    0
#>          features
#> docs      forcing thousands
#>   art_001       1         1
#>   art_002       0         0
#>   art_003       0         1
#>   art_004       0         0
#>   art_005       0         0
#>   art_006       0         0
#> [ reached max_ndoc ... 14 more documents, reached max_nfeat ... 243 more features ]

set.seed(42)
textplot_wordcloud(dfm_articles, max_words = 80,
                   color = c("#1e3a2f", "#4a7c5e", "#8aad96", 
                             "#c8832a", "#7a3a2f"))

Word cloud of the simulated corpus. Larger words appear more frequently across all 20 articles. Environmental terms (wildfire, water, community, climate) dominate.

Available dictionaries in R

Several dictionaries are available through different R packages. Each was built for a different purpose.

Scoring with LSD2015 (recommended for news)

The Lexicoder Sentiment Dictionary (LSD2015) was designed specifically for political and news text. It includes negation-aware word patterns (e.g., “not good” is captured as negative). It is built into quanteda as data_dictionary_LSD2015.

# Apply LSD2015 dictionary to the dfm
dfm_sent <- dfm_lookup(dfm_articles,
                        dictionary = data_dictionary_LSD2015)

# Convert to a tidy data frame
lsd_scores <- convert(dfm_sent, to = "data.frame") %>%
  rename(doc_id = doc_id) %>%
  mutate(
    net_sentiment = positive - negative,
    total_words   = ntoken(corp),
    sent_ratio    = net_sentiment / total_words
  ) %>%
  bind_cols(docvars(corp))

head(lsd_scores)

#>    doc_id negative positive neg_positive neg_negative net_sentiment total_words
#> 1 art_001        3        1            0            0            -2          20
#> 2 art_002        3        3            0            0             0          24
#> 3 art_003        0        2            0            0             2          26
#> 4 art_004        3        0            0            0            -3          18
#> 5 art_005        1        0            0            0            -1          22
#> 6 art_006        3        1            0            0            -2          19
#>    sent_ratio           source       date        topic
#> 1 -0.10000000 National Tribune 2024-01-15     wildfire
#> 2  0.00000000 National Tribune 2024-02-15        water
#> 3  0.07692308 National Tribune 2024-03-15       energy
#> 4 -0.16666667 National Tribune 2024-04-15 biodiversity
#> 5 -0.04545455 National Tribune 2024-05-15      climate
#> 6 -0.10526316   Valley Gazette 2024-06-15     wildfire

ggplot(lsd_scores, aes(x = source, y = net_sentiment, fill = source)) +
  geom_boxplot(width = 0.45, outlier.shape = NA,
               colour = "#1e3a2f", alpha = 0.75, linewidth = 0.6) +
  geom_jitter(aes(colour = source), width = 0.15,
              size = 2, alpha = 0.6) +
  geom_hline(yintercept = 0, linetype = "dashed", colour = "#7a7a72") +
  scale_fill_manual(values = c(
    "Coastal Herald"     = "#4a7c5e",
    "Mountain Post"      = "#7a3a2f",
    "National Tribune"   = "#1e3a2f",
    "Valley Gazette"     = "#c8832a"
  )) +
  scale_colour_manual(values = c(
    "Coastal Herald"     = "#4a7c5e",
    "Mountain Post"      = "#7a3a2f",
    "National Tribune"   = "#1e3a2f",
    "Valley Gazette"     = "#c8832a"
  )) +
  labs(x = NULL, y = "Net sentiment (positive - negative words)",
       title = "LSD2015 Sentiment by News Source") +
  theme_minimal(base_size = 12) +
  theme(
    plot.title         = element_text(family = "serif", face = "bold",
                                      colour = "#1e3a2f"),
    legend.position    = "none",
    panel.grid.major.x = element_blank()
  )

Net LSD2015 sentiment (positive minus negative word counts) by news source. Each point is one article; the dashed line at zero separates net-positive from net-negative.

Scoring with tidytext and multiple lexicons

The tidytext package provides access to AFINN, Bing, and NRC lexicons through a tidy one-word-per-row workflow. Running multiple lexicons in parallel is the best way to check whether your findings are robust to dictionary choice.

# Tidy tokenization
tidy_articles <- articles %>%
  unnest_tokens(word, text)

# Score with AFINN
afinn_scores <- tidy_articles %>%
  inner_join(get_sentiments("afinn"), by = "word") %>%
  group_by(doc_id, source, date, topic) %>%
  summarize(
    afinn_mean = mean(value),
    afinn_sum  = sum(value),
    n_scored   = n(),
    .groups = "drop"
  )

# Score with Bing
bing_scores <- tidy_articles %>%
  inner_join(get_sentiments("bing"), by = "word") %>%
  count(doc_id, source, date, topic, sentiment) %>%
  pivot_wider(names_from = sentiment, values_from = n, values_fill = 0) %>%
  mutate(bing_net = positive - negative)

# Score with NRC (emotion categories)
nrc_emotions <- tidy_articles %>%
  inner_join(get_sentiments("nrc"), by = "word") %>%
  count(doc_id, source, topic, sentiment) %>%
  pivot_wider(names_from = sentiment, values_from = n, values_fill = 0)

nrc_topic <- tidy_articles %>%
  inner_join(get_sentiments("nrc"), by = "word") %>%
  filter(!sentiment %in% c("positive", "negative")) %>%
  count(topic, sentiment)

ggplot(nrc_topic, aes(x = reorder(sentiment, n), y = n, fill = topic)) +
  geom_col(position = "dodge", alpha = 0.85) +
  coord_flip() +
  scale_fill_manual(values = c(
    "biodiversity" = "#4a7c5e",
    "climate"      = "#7a3a2f",
    "energy"       = "#c8832a",
    "water"        = "#2d5a3d",
    "wildfire"     = "#8b3a2a"
  )) +
  labs(x = NULL, y = "Word count",
       title = "NRC Emotion Profiles by Topic",
       fill = "Topic") +
  theme_minimal(base_size = 12) +
  theme(
    plot.title = element_text(family = "serif", face = "bold",
                               colour = "#1e3a2f"),
    legend.position = "bottom"
  )

NRC emotion profile aggregated by topic across all articles. Each bar represents the total word count matching that emotion category. Fear and negative words dominate wildfire and climate articles; trust and positive words are more common in energy and biodiversity articles.

Handling valence shifters with sentimentr

A critical weakness of bag-of-words lexicons is that they score “not good” as positive (the word “good” matches the dictionary). The sentimentr package handles valence shifters – negators (“not”), amplifiers (“very”), de-amplifiers (“somewhat”), and adversative conjunctions (“but”) – at the sentence level.

# sentimentr works on raw text, not tokenized input
sr_scores <- sentiment_by(articles$text, by = articles$doc_id)

# Merge back with metadata
sr_df <- articles %>%
  left_join(sr_scores, by = c("doc_id" = "doc_id"))

head(sr_df %>% select(doc_id, source, topic, ave_sentiment))

#> # A tibble: 6 × 4
#>   doc_id  source           topic        ave_sentiment
#>   <chr>   <chr>            <chr>                <dbl>
#> 1 art_001 National Tribune wildfire           -0.327 
#> 2 art_002 National Tribune water               0.341 
#> 3 art_003 National Tribune energy              0.338 
#> 4 art_004 National Tribune biodiversity       -0.449 
#> 5 art_005 National Tribune climate             0.0671
#> 6 art_006 Valley Gazette   wildfire           -0.342

# Sentence-level detail for two articles
detail <- sentiment(articles$text[c(1, 3)])
detail$article <- rep(c("art_001 (Wildfire)", "art_003 (Solar)"),
                      times = c(sum(detail$element_id == 1),
                                sum(detail$element_id == 2)))

ggplot(detail, aes(x = sentence_id, y = sentiment, fill = sentiment > 0)) +
  geom_col(alpha = 0.8, show.legend = FALSE) +
  facet_wrap(~ article, scales = "free_x") +
  scale_fill_manual(values = c("TRUE" = "#4a7c5e", "FALSE" = "#7a3a2f")) +
  labs(x = "Sentence", y = "Sentiment score",
       title = "Sentence-Level Sentiment (sentimentr)") +
  theme_minimal(base_size = 11) +
  theme(
    plot.title = element_text(family = "serif", face = "bold",
                               colour = "#1e3a2f"),
    strip.text = element_text(face = "bold")
  )

Sentence-level sentiment from sentimentr for two example articles, one about wildfire (article 1) and one about solar energy (article 3). The valence-shifted scoring correctly identifies the overall negative tone of wildfire coverage and positive tone of the solar story.

Targeted dictionary lookup in quanteda

The simplest ABSA approach uses quanteda’s windowed token selection. You define your target entity, extract a window of +/- 10 tokens around every mention of that entity, and then apply a sentiment dictionary only to the text within that window.

# Define target entities
targets <- list(
  wildfire  = c("wildfire", "fire", "burn", "blaze"),
  renewable = c("solar", "wind", "clean energy", "renewable"),
  water     = c("water", "reservoir", "drought", "runoff")
)

# Function: extract sentiment around a target
target_sentiment <- function(toks, target_words, window = 10) {
  # Select tokens within +/- window of target
  toks_target <- tokens_keep(toks, pattern = target_words,
                              window = window, padding = FALSE)
  
  # Apply LSD2015
  dfm_target <- dfm(toks_target) %>%
    dfm_lookup(dictionary = data_dictionary_LSD2015)
  
  result <- convert(dfm_target, to = "data.frame") %>%
    mutate(
      net = positive - negative,
      target = paste(target_words, collapse = "/")
    )
  
  result
}

# Score each target
absa_results <- map_dfr(targets, ~ target_sentiment(toks, .x)) %>%
  bind_cols(
    docvars(corp) %>% slice(rep(1:n(), times = length(targets)))
  )

# Summarize by target
absa_summary <- absa_results %>%
  group_by(target) %>%
  summarize(
    mean_net   = mean(net, na.rm = TRUE),
    median_net = median(net, na.rm = TRUE),
    sd_net     = sd(net, na.rm = TRUE),
    n_docs     = n(),
    .groups = "drop"
  )

kable(absa_summary, digits = 2,
      col.names = c("Target Entity", "Mean Net", "Median Net", "SD", "N"),
      caption = "Aspect-level sentiment: average net LSD2015 sentiment in a +/- 10 word window around each target entity.") %>%
  kable_styling(full_width = FALSE, font_size = 14)

ClimateBERT: a domain-adapted model

ClimateBERT is a DistilRoBERTa model further pre-trained on over 2 million climate-related paragraphs from news, research articles, and corporate reports. The sentiment variant (climatebert/distilroberta-base-climate-sentiment) classifies text into three categories: opportunity, neutral, and risk. On its benchmark dataset, it reports weighted F1 of approximately 0.84.

Calling transformers from R with the text package

The text package provides an R bridge to HuggingFace models. It requires a one-time Python environment setup (handled by textrpp_install()).

# ── One-time setup (run once, not every session) ───────────────────────
# install.packages("text")
library(text)
textrpp_install()        # sets up a conda environment with PyTorch
textrpp_initialize()     # connects R to the Python backend

# ── Classify articles with ClimateBERT ─────────────────────────────────
library(text)

# Use the climate-sentiment model
climate_sent <- textClassify(
  articles$text,
  model = "climatebert/distilroberta-base-climate-sentiment",
  device = "cpu",             # use "cuda" if GPU available
  return_incorrect_results = FALSE
)

# Result: a data frame with columns for each class probability
# and a predicted label
head(climate_sent)

Other models worth knowing

When transformers are not worth the cost

If your hand-coded validation shows that your lexicon approach achieves macro-F1 above 0.65 and human inter-annotator agreement (Cohen’s kappa) is below 0.55, adding a transformer layer will not produce meaningful gains. The annotation noise itself becomes the ceiling. Invest the time in better annotation guidelines instead.

The three aggregation dimensions

Building a sentiment time series involves three sequential choices:

1. Within-document: How do you combine sentence-level scores into one document score? Options include the mean, median, proportional weighting (longer sentences count more), or max-absolute (capture peak intensity).

2. Across-documents per period: How do you combine document scores within a week or month? Equal weighting treats every article the same. Article-count proportional weighting lets outlets with more articles have more influence. Source-balanced weighting averages within each outlet first, then across outlets.

3. Across-time smoothing: Raw weekly scores are noisy. Smoothing options include moving averages, exponential decay, and Beta-weighted curves that emphasize recent observations.

Building a simple monthly index

# Monthly sentiment index from LSD2015 scores
monthly <- lsd_scores %>%
  mutate(month = floor_date(date, "month")) %>%
  group_by(month) %>%
  summarize(
    mean_sentiment = mean(sent_ratio, na.rm = TRUE),
    sd_sentiment   = sd(sent_ratio, na.rm = TRUE),
    n_articles     = n(),
    se             = sd_sentiment / sqrt(n_articles),
    lower          = mean_sentiment - 1.96 * se,
    upper          = mean_sentiment + 1.96 * se,
    .groups = "drop"
  )

ggplot(monthly, aes(x = month, y = mean_sentiment)) +
  geom_ribbon(aes(ymin = lower, ymax = upper),
              fill = "#8aad96", alpha = 0.3) +
  geom_line(colour = "#1e3a2f", linewidth = 0.9) +
  geom_point(colour = "#1e3a2f", size = 2.5) +
  geom_hline(yintercept = 0, linetype = "dashed", colour = "#7a7a72") +
  scale_x_date(date_labels = "%b %Y") +
  labs(x = NULL, y = "Mean sentiment ratio",
       title = "Monthly Environmental News Sentiment Index",
       subtitle = "LSD2015 dictionary, 95% CI ribbon") +
  theme_minimal(base_size = 12) +
  theme(
    plot.title    = element_text(family = "serif", face = "bold",
                                 colour = "#1e3a2f"),
    plot.subtitle = element_text(colour = "#7a7a72", size = 10)
  )

Monthly sentiment index from LSD2015 dictionary scoring. Ribbon shows 95% confidence intervals. With only 20 simulated articles, the intervals are wide; a real project would have hundreds or thousands of articles per month.

The sentometrics package

For a more principled approach, the sentometrics package (Ardia et al., Journal of Statistical Software 2021) handles all three aggregation dimensions in a single workflow. It also provides elastic-net regression tools for relating sentiment indices to external variables.

# ── sentometrics workflow (conceptual) ─────────────────────────────────
library(sentometrics)

# 1. Create a sentocorpus (requires 'id', 'date', 'texts' columns,
#    plus optional feature columns for topic or source)
sc <- sento_corpus(articles,
                   corpusdf = TRUE,
                   do.clean = TRUE)

# 2. Define lexicons
lexicons <- sento_lexicons(list(
  LSD2015 = data_dictionary_LSD2015,   # from quanteda
  AFINN   = get_sentiments("afinn")     # from tidytext
))

# 3. Define aggregation controls
ctr <- ctr_agg(
  howWithin  = "proportional",     # weight by sentence length
  howDocs    = "equal_weight",     # equal across documents
  howTime    = "equal_weight",     # equal across time
  by         = "month",
  lag        = 6                   # 6-month rolling window
)

# 4. Compute sentiment measures
sento_measures <- sento_measures(sc, lexicons, ctr)

# 5. Plot
plot(sento_measures)

Step 1: Hand-code a sample

Draw a stratified random sample of 200 to 500 articles (stratified by source, time period, and topic). Have two independent coders classify each article on the same scale your model uses (e.g., positive / neutral / negative, or opportunity / neutral / risk).

Step 2: Measure human agreement

Before comparing the model to humans, measure how well the humans agree with each other. If the coders disagree substantially, the task itself is too subjective for any model to solve reliably.

# Simulated coding data: two coders, 50 articles, 3-class scale
set.seed(123)
labels <- c("positive", "neutral", "negative")
coder1 <- sample(labels, 50, replace = TRUE, prob = c(0.25, 0.45, 0.30))
coder2 <- coder1   # start with agreement, then add noise
noise_idx <- sample(50, 12)
coder2[noise_idx] <- sample(labels, 12, replace = TRUE)

# Cohen's kappa
confusion <- table(Coder1 = coder1, Coder2 = coder2)
confusion

#>           Coder2
#> Coder1     negative neutral positive
#>   negative       11       2        2
#>   neutral         1      19        1
#>   positive        2       0       12

# Manual Cohen's kappa calculation
n <- sum(confusion)
p_o <- sum(diag(confusion)) / n               # observed agreement
p_e <- sum(rowSums(confusion) * colSums(confusion)) / n^2  # expected by chance
kappa <- (p_o - p_e) / (1 - p_e)

cat("Observed agreement (p_o):", round(p_o, 3), "\n")

#> Observed agreement (p_o): 0.84

cat("Expected agreement (p_e):", round(p_e, 3), "\n")

#> Expected agreement (p_e): 0.344

cat("Cohen's kappa:           ", round(kappa, 3), "\n")

#> Cohen's kappa:            0.756

Step 3: Compare model to human consensus

After the two coders resolve disagreements (creating a consensus label), compute the same metrics between the model’s predictions and the consensus labels.

# Simulated model predictions vs consensus labels
consensus <- coder1  # in practice, resolved disagreements
model_pred <- consensus
# Introduce some model errors
model_err <- sample(50, 8)
model_pred[model_err] <- sample(labels, 8, replace = TRUE)

# Confusion matrix
conf_model <- table(Model = model_pred, Human = consensus)

# Per-class precision, recall, F1
calc_metrics <- function(conf_mat) {
  classes <- rownames(conf_mat)
  metrics <- map_dfr(classes, function(cls) {
    tp <- conf_mat[cls, cls]
    fp <- sum(conf_mat[cls, ]) - tp
    fn <- sum(conf_mat[, cls]) - tp
    precision <- ifelse(tp + fp > 0, tp / (tp + fp), 0)
    recall    <- ifelse(tp + fn > 0, tp / (tp + fn), 0)
    f1        <- ifelse(precision + recall > 0,
                        2 * precision * recall / (precision + recall), 0)
    tibble(Class = cls, Precision = precision,
           Recall = recall, F1 = f1)
  })
  metrics
}

per_class <- calc_metrics(conf_model)
macro_f1 <- mean(per_class$F1)

kable(per_class, digits = 3,
      caption = paste0("Per-class metrics (model vs. human consensus). ",
                       "Macro-F1 = ", round(macro_f1, 3), ".")) %>%
  kable_styling(full_width = FALSE, font_size = 14)

What metrics to report

Anchoring sentiment to environmental events

One of the most informative analyses is to overlay your sentiment time series with known environmental events. This reveals whether sentiment shifts lead, lag, or coincide with real-world triggers.

# Simulated event markers
events <- tibble(
  date  = as.Date(c("2024-03-15", "2024-07-15", "2024-11-15")),
  event = c("IPCC Report\nRelease", "Record Wildfire\nSeason Begins", 
            "COP29\nSummit")
)

ggplot(monthly, aes(x = month, y = mean_sentiment)) +
  geom_ribbon(aes(ymin = lower, ymax = upper),
              fill = "#8aad96", alpha = 0.25) +
  geom_line(colour = "#1e3a2f", linewidth = 0.9) +
  geom_point(colour = "#1e3a2f", size = 2.5) +
  geom_vline(data = events, aes(xintercept = date),
             linetype = "dotted", colour = "#c8832a", linewidth = 0.7) +
  geom_label(data = events, aes(x = date, y = max(monthly$upper, na.rm = TRUE),
                                 label = event),
             size = 2.5, fill = "#fdf8f0", colour = "#c8832a",
             label.padding = unit(0.3, "lines"), family = "mono") +
  geom_hline(yintercept = 0, linetype = "dashed", colour = "#7a7a72") +
  scale_x_date(date_labels = "%b %Y") +
  labs(x = NULL, y = "Mean sentiment ratio",
       title = "Environmental News Sentiment with Event Markers") +
  theme_minimal(base_size = 12) +
  theme(
    plot.title = element_text(family = "serif", face = "bold",
                               colour = "#1e3a2f")
  )

Sentiment time series with event annotations. In a real analysis with hundreds of articles, you would mark COP summits, IPCC report releases, major disasters, and policy announcements.

The three approaches compared

We score every article with three different methods and compare their behavior. This is exactly the multi-lexicon robustness check recommended in Section 4, but now we have an external criterion (price movements) to evaluate against.

Building an energy news corpus

In a production workflow, you would accumulate articles over weeks or months from RSS feeds (Reuters Energy, CNBC Energy, OilPrice.com, EIA Today in Energy) and store them in a cache file. For this tutorial, we build a synthetic corpus of representative energy headlines that covers the same topic channels the real feeds would provide.

# ── Synthetic energy news corpus ───────────────────────────────────────
# In practice, these would come from RSS feeds via xml2::read_xml()
# or from a cached CSV built by a scheduled scraping script.
set.seed(42)

headlines_pool <- c(
  # Supply risk (negative)
  "Oil prices surge as Middle East tensions escalate after military strikes",
  "Pipeline attack in Nigeria disrupts 200,000 barrels per day of supply",
  "Conflict in Libya shuts major oil field cutting output by 300,000 bpd",
  "Geopolitical tensions push crude above $100 as supply risks mount",
  "War risk premium builds as military forces mobilize near oil infrastructure",
  "Embargo threat from major producers sends Brent to six-month high",
  # Supply risk (positive)
  "Ceasefire agreement reached easing supply disruption fears in Middle East",
  "Pipeline repairs completed restoring full capacity after weeks of disruption",
  # OPEC policy
  "OPEC+ agrees to extend production cuts through end of year",
  "Saudi Arabia announces voluntary output reduction of 1 million barrels per day",
  "OPEC meeting ends with surprise decision to maintain current quotas",
  "UAE pushes for higher production target at contentious OPEC+ session",
  "OPEC+ considers easing output curbs as oil prices stabilize above $80",
  # Demand / macroeconomic

  "Global oil demand forecast revised upward on strong Chinese economic recovery",
  "US GDP growth slows sharply raising recession fears and oil demand concerns",
  "Strong travel demand pushes jet fuel consumption above pre-pandemic levels",
  "Economic recession fears mount as consumer spending falls sharply",
  "Driving season kicks off with record US gasoline consumption forecast",
  "Industrial output in Europe falls for second quarter amid energy cost pressures",
  # Refinery
  "Major refinery fire in Texas shuts 250,000 bpd capacity for weeks",
  "Spring maintenance season begins with multiple refineries taking units offline",
  "Heating oil stocks at five-year low ahead of winter demand season",
  "Crack spread widens as refinery capacity offline exceeds seasonal norms",
  # Mixed / neutral
  "Oil markets steady as traders weigh competing supply and demand signals",
  "Energy stocks mixed as crude prices trade in narrow range near $75 per barrel",
  "US crude inventories rise more than expected according to EIA weekly report",
  "Analysts divided on crude oil direction as technical signals conflict"
)

# Generate 300 articles across 18 months
DATE_START <- Sys.Date() - 540
DATE_END   <- Sys.Date()
n_articles <- 300

energy_news <- tibble(
  date   = sample(seq.Date(DATE_START + 7, DATE_END - 7, by = "day"),
                  n_articles, replace = TRUE) |> sort(),
  text   = sample(headlines_pool, n_articles, replace = TRUE),
  source = sample(c("Reuters_Energy", "CNBC_Energy", "EIA_Today"),
                  n_articles, replace = TRUE)
)

cat(sprintf("Energy corpus: %d articles (%s to %s)\n",
            nrow(energy_news), min(energy_news$date), max(energy_news$date)))

#> Energy corpus: 300 articles (2024-11-24 to 2026-05-03)

Scoring with three approaches simultaneously

This is the core comparison. We apply all three methods to the same articles, producing three parallel sentiment columns.

# ── Approach 1: Syuzhet (Stanford NLP lexicon) ─────────────────────────
# Returns a numeric score per text element. No negation handling.
energy_news$sent_syuzhet <- syuzhet::get_sentiment(
  energy_news$text, method = "syuzhet"
)

# ── Approach 2: LSD2015 via quanteda ──────────────────────────────────
# Purpose-built for news. Includes negation patterns.
# Uses tokens_lookup() to capture multi-word patterns before dfm.
energy_corp <- corpus(energy_news$text)
energy_toks <- tokens(energy_corp, remove_punct = TRUE,
                      remove_numbers = TRUE) |>
  tokens_tolower()

toks_lsd <- tokens_lookup(energy_toks,
                          dictionary = data_dictionary_LSD2015,
                          exclusive = TRUE)
lsd_dfm  <- dfm(toks_lsd)
lsd_mat  <- as.matrix(lsd_dfm)
total_toks <- ntoken(energy_toks)

# The LSD2015 has four keys: negative, positive, neg_negative, neg_positive
safe_col <- function(mat, col_name) {
  if (col_name %in% colnames(mat)) mat[, col_name] else rep(0, nrow(mat))
}
pos_signal <- safe_col(lsd_mat, "positive") + safe_col(lsd_mat, "neg_negative")
neg_signal <- safe_col(lsd_mat, "negative") + safe_col(lsd_mat, "neg_positive")

energy_news$sent_lsd2015 <- (pos_signal - neg_signal) / pmax(total_toks, 1)

# ── Approach 3: sentimentr (valence-shifted) ──────────────────────────
# Handles negation, amplifiers, and adversative conjunctions.
clean_text <- energy_news$text |>
  gsub("<[^>]+>", " ", x = _) |>
  gsub("[[:cntrl:]]", " ", x = _) |>
  gsub("\\s+", " ", x = _) |>
  trimws()
clean_text <- ifelse(grepl("[.!?]\\s*$", clean_text),
                     clean_text,
                     paste0(clean_text, "."))

sr_scores <- sentiment_by(get_sentences(clean_text))
energy_news$sent_sentimentr <- sr_scores$ave_sentiment

How well do the three approaches agree?

Before looking at price correlations, check whether the methods even agree with each other. If they produce wildly different rankings, that tells you the “sentiment” construct is not stable across measurement instruments.

# Inter-approach correlations
approach_cor <- cor(
  energy_news |> select(sent_syuzhet, sent_lsd2015, sent_sentimentr),
  use = "pairwise.complete.obs"
)

kable(round(approach_cor, 3),
      caption = "Pairwise Pearson correlations between three sentiment approaches on the same energy news corpus.") %>%
  kable_styling(full_width = FALSE, font_size = 14)

energy_news |>
  select(sent_syuzhet, sent_lsd2015, sent_sentimentr) |>
  pivot_longer(everything(), names_to = "approach", values_to = "score") |>
  mutate(approach = recode(approach,
    "sent_syuzhet"   = "Syuzhet",
    "sent_lsd2015"   = "LSD2015",
    "sent_sentimentr" = "sentimentr"
  )) |>
  ggplot(aes(approach, score, fill = approach)) +
  geom_boxplot(alpha = 0.75, outlier.size = 0.8, width = 0.5,
               colour = "#1e3a2f") +
  geom_hline(yintercept = 0, linewidth = 0.3, colour = "#7a7a72") +
  scale_fill_manual(values = c(
    "Syuzhet"    = "#2166AC",
    "LSD2015"    = "#B2182B",
    "sentimentr" = "#1B7837"
  )) +
  labs(x = NULL, y = "Sentiment score",
       title = "Score Distributions by Approach") +
  theme_minimal(base_size = 12) +
  theme(
    plot.title      = element_text(family = "serif", face = "bold",
                                    colour = "#1e3a2f"),
    legend.position = "none"
  )

Score distributions by approach. Syuzhet produces integer-valued scores with wider spread. LSD2015 is proportion-based and concentrated near zero. sentimentr falls between the two.

Topic-segmented sentiment

Energy news is not monolithic. An article about OPEC production cuts has a different sentiment profile than an article about a refinery fire, even though both affect prices. We classify articles into four topic channels using keyword matching and score each channel separately.

# Topic keyword dictionaries
topic_keywords <- list(
  supply_risk = c("sanctions", "geopolit", "middle east", "pipeline attack",
                  "conflict", "war", "embargo", "disruption", "ceasefire"),
  opec_policy = c("opec", "production cut", "output", "quota",
                  "saudi", "uae", "barrel per day"),
  demand      = c("recession", "gdp", "demand", "economic",
                  "consumption", "driving season", "industrial output"),
  refinery    = c("refinery", "outage", "shutdown", "maintenance",
                  "heating oil", "crack spread", "capacity")
)

# Tag each article
text_lower <- tolower(energy_news$text)
for (topic in names(topic_keywords)) {
  pattern <- paste(topic_keywords[[topic]], collapse = "|")
  energy_news[[paste0("is_", topic)]] <- grepl(pattern, text_lower)
}

# Topic hit summary
topic_counts <- energy_news |>
  summarize(across(starts_with("is_"), sum)) |>
  pivot_longer(everything(), names_to = "topic", values_to = "n") |>
  mutate(topic = gsub("^is_", "", topic),
         pct = round(100 * n / nrow(energy_news), 1))

kable(topic_counts, col.names = c("Topic", "Articles", "% of Corpus"),
      caption = "Topic channel coverage in the energy news corpus.") %>%
  kable_styling(full_width = FALSE, font_size = 14)

topic_sent <- energy_news |>
  pivot_longer(starts_with("is_"), names_to = "topic", values_to = "hit") |>
  filter(hit) |>
  mutate(topic = gsub("^is_", "", topic) |>
           gsub("_", " ", x = _) |> tools::toTitleCase()) |>
  group_by(topic) |>
  summarize(
    mean_sent = mean(sent_syuzhet),
    se        = sd(sent_syuzhet) / sqrt(n()),
    n         = n(),
    .groups   = "drop"
  )

ggplot(topic_sent, aes(reorder(topic, mean_sent), mean_sent, fill = mean_sent > 0)) +
  geom_col(alpha = 0.8, width = 0.6, colour = "#1e3a2f", show.legend = FALSE) +
  geom_errorbar(aes(ymin = mean_sent - se, ymax = mean_sent + se),
                width = 0.15, linewidth = 0.7, colour = "#1e3a2f") +
  geom_hline(yintercept = 0, linewidth = 0.4) +
  scale_fill_manual(values = c("TRUE" = "#4a7c5e", "FALSE" = "#7a3a2f")) +
  coord_flip() +
  labs(x = NULL, y = "Mean Syuzhet sentiment",
       title = "Sentiment by Topic Channel") +
  theme_minimal(base_size = 12) +
  theme(
    plot.title = element_text(family = "serif", face = "bold",
                               colour = "#1e3a2f")
  )

Mean Syuzhet sentiment by topic channel. Supply risk articles are the most negative; articles about demand recovery tend to be mixed. Error bars show +/- 1 standard error.

Aggregating to weekly and merging with price data

To correlate sentiment with commodity prices, we need both series on the same time scale. EIA petroleum prices are published weekly (Wednesday, for the prior week), so we aggregate daily sentiment to ISO weeks.

# Weekly sentiment aggregation
energy_weekly <- energy_news |>
  mutate(week_date = floor_date(date, "week", week_start = 1)) |>
  group_by(week_date) |>
  summarize(
    syuzhet_mean    = mean(sent_syuzhet, na.rm = TRUE),
    lsd2015_mean    = mean(sent_lsd2015, na.rm = TRUE),
    sentimentr_mean = mean(sent_sentimentr, na.rm = TRUE),
    n_articles      = n(),
    pct_negative    = mean(sent_syuzhet < 0, na.rm = TRUE),
    .groups = "drop"
  )

cat(sprintf("Weekly aggregation: %d weeks, %.1f articles/week on average\n",
            nrow(energy_weekly), mean(energy_weekly$n_articles)))

#> Weekly aggregation: 74 weeks, 4.1 articles/week on average

# ── Simulated weekly WTI crude prices ──────────────────────────────────
# In the full workflow, these come from fetch_eia_v2() in the companion script.
# We simulate a plausible price series with trend + noise.
set.seed(99)
n_weeks <- nrow(energy_weekly)
trend <- seq(72, 85, length.out = n_weeks) +
  cumsum(rnorm(n_weeks, 0, 1.5))

price_weekly <- tibble(
  week_date = energy_weekly$week_date,
  wti_spot  = round(pmax(trend, 55), 2)
) |>
  mutate(wti_return = (wti_spot / lag(wti_spot) - 1))

# Merge sentiment + prices
merged_energy <- left_join(energy_weekly, price_weekly, by = "week_date") |>
  drop_na(wti_return)

Correlation analysis: which approach tracks prices best?

Now we have three sentiment columns and one price-return column on the same weekly grid. The question is: which sentiment approach has the strongest contemporaneous correlation with price movements?

# Contemporaneous correlations
approach_cols <- c("syuzhet_mean", "lsd2015_mean", "sentimentr_mean")

cor_results <- tibble(approach = approach_cols) |>
  rowwise() |>
  mutate(
    r       = cor(merged_energy[[approach]], merged_energy$wti_return,
                  use = "pairwise.complete.obs"),
    p_value = cor.test(merged_energy[[approach]],
                       merged_energy$wti_return)$p.value,
    label   = recode(approach,
      "syuzhet_mean"    = "Syuzhet",
      "lsd2015_mean"    = "LSD2015",
      "sentimentr_mean" = "sentimentr"
    )
  ) |>
  ungroup() |>
  mutate(
    sig = case_when(
      p_value < 0.01 ~ "***",
      p_value < 0.05 ~ "**",
      p_value < 0.10 ~ "*",
      TRUE ~ ""
    )
  )

kable(cor_results |> select(label, r, p_value, sig),
      digits = 4,
      col.names = c("Approach", "r", "p-value", "Sig."),
      caption = "Contemporaneous Pearson correlations: weekly mean sentiment vs. weekly WTI return. On synthetic data, these correlations are near zero (as expected). With real news and real prices, modest correlations (|r| = 0.05-0.20) are typical.") %>%
  kable_styling(full_width = FALSE, font_size = 14)

Lead-lag analysis: does sentiment predict price?

Contemporaneous correlation does not tell you about predictive value. A cross-correlation function (CCF) tests whether sentiment at time \(t - k\) correlates with price returns at time \(t\). A negative lag with a strong CCF means sentiment leads price, suggesting some predictive content.

# CCF for Syuzhet vs WTI returns
valid <- complete.cases(merged_energy$syuzhet_mean, merged_energy$wti_return)
if (sum(valid) >= 20) {
  cc <- ccf(merged_energy$syuzhet_mean[valid],
            merged_energy$wti_return[valid],
            lag.max = 4, plot = FALSE)
  
  ccf_df <- tibble(
    lag = cc$lag[, 1, 1],
    ccf = cc$acf[, 1, 1]
  )
  
  ggplot(ccf_df, aes(lag, ccf)) +
    geom_hline(yintercept = 0, linewidth = 0.3) +
    geom_vline(xintercept = 0, linetype = "dashed", colour = "#7a7a72") +
    geom_col(aes(fill = lag < 0), alpha = 0.75, width = 0.6,
             show.legend = FALSE) +
    scale_fill_manual(values = c("TRUE" = "#2166AC", "FALSE" = "#B2182B")) +
    labs(x = "Lag (weeks, negative = sentiment leads)",
         y = "Cross-correlation",
         title = "Lead-Lag: Syuzhet Sentiment vs. WTI Returns") +
    theme_minimal(base_size = 12) +
    theme(
      plot.title = element_text(family = "serif", face = "bold",
                                 colour = "#1e3a2f")
    )
}

Cross-correlation between Syuzhet sentiment and WTI weekly returns at lags of -4 to +4 weeks. Negative lags (left side) represent sentiment leading price. With synthetic data the signal is weak; real data from the companion script typically shows modest leading correlations.

Diagnostic multi-panel figure

The companion script produces an 8-panel diagnostic figure. Here we build a simpler 4-panel version that captures the key comparisons.

approach_colors <- c(
  "Syuzhet"    = "#2166AC",
  "LSD2015"    = "#B2182B",
  "sentimentr" = "#1B7837"
)

# Panel 1: Sentiment time series
p1_data <- merged_energy |>
  select(week_date, syuzhet_mean, lsd2015_mean, sentimentr_mean) |>
  pivot_longer(-week_date, names_to = "approach", values_to = "score") |>
  mutate(approach = recode(approach,
    "syuzhet_mean"    = "Syuzhet",
    "lsd2015_mean"    = "LSD2015",
    "sentimentr_mean" = "sentimentr"
  ))

p1 <- ggplot(p1_data, aes(week_date, score, colour = approach)) +
  geom_hline(yintercept = 0, linewidth = 0.3, colour = "#7a7a72") +
  geom_line(linewidth = 0.4, alpha = 0.5) +
  geom_smooth(method = "loess", se = FALSE, linewidth = 0.9, span = 0.3) +
  scale_colour_manual(values = approach_colors) +
  labs(title = "Weekly Sentiment: Three Approaches",
       x = NULL, y = "Mean score", colour = NULL) +
  theme_minimal(base_size = 10) +
  theme(plot.title      = element_text(family = "serif", face = "bold",
                                        colour = "#1e3a2f"),
        legend.position = "bottom")

# Panel 2: WTI with sentiment overlay (scaled)
wti_range  <- range(merged_energy$wti_spot, na.rm = TRUE)
sent_range <- range(merged_energy$syuzhet_mean, na.rm = TRUE)
scale_fn <- function(s) {
  (s - sent_range[1]) / max(diff(sent_range), 0.001) *
    diff(wti_range) * 0.3 + wti_range[1] - diff(wti_range) * 0.05
}

p2 <- ggplot(merged_energy, aes(week_date)) +
  geom_line(aes(y = wti_spot), colour = "black", linewidth = 0.7) +
  geom_line(aes(y = scale_fn(syuzhet_mean)),
            colour = "#2166AC", alpha = 0.5) +
  geom_smooth(aes(y = scale_fn(syuzhet_mean)),
              method = "loess", se = FALSE, span = 0.3,
              colour = "#2166AC", linewidth = 0.8) +
  labs(title = "WTI Crude vs. Syuzhet Sentiment (scaled)",
       subtitle = "Black = WTI ($/bbl), blue = sentiment (scaled)",
       x = NULL, y = "WTI Spot ($/bbl)") +
  theme_minimal(base_size = 10) +
  theme(plot.title    = element_text(family = "serif", face = "bold",
                                      colour = "#1e3a2f"),
        plot.subtitle = element_text(size = 8, colour = "#7a7a72"))

# Panel 3: Distributions (reuse from earlier)
p3_data <- energy_news |>
  select(sent_syuzhet, sent_lsd2015, sent_sentimentr) |>
  pivot_longer(everything(), names_to = "approach", values_to = "score") |>
  mutate(approach = recode(approach,
    "sent_syuzhet"    = "Syuzhet",
    "sent_lsd2015"    = "LSD2015",
    "sent_sentimentr" = "sentimentr"
  ))

p3 <- ggplot(p3_data, aes(approach, score, fill = approach)) +
  geom_boxplot(alpha = 0.7, outlier.size = 0.6, width = 0.5,
               colour = "#1e3a2f") +
  geom_hline(yintercept = 0, linewidth = 0.3) +
  scale_fill_manual(values = approach_colors) +
  labs(title = "Score Distributions", x = NULL, y = "Score") +
  theme_minimal(base_size = 10) +
  theme(plot.title      = element_text(family = "serif", face = "bold",
                                        colour = "#1e3a2f"),
        legend.position = "none")

# Panel 4: Article volume
p4 <- ggplot(merged_energy, aes(week_date, n_articles)) +
  geom_col(fill = "#8aad96", alpha = 0.7) +
  geom_smooth(method = "loess", se = FALSE, colour = "#1e3a2f",
              span = 0.3, linewidth = 0.7) +
  labs(title = "Weekly Article Volume",
       x = NULL, y = "Articles / week") +
  theme_minimal(base_size = 10) +
  theme(plot.title = element_text(family = "serif", face = "bold",
                                    colour = "#1e3a2f"))

# Assemble with patchwork
(p1 + p2) / (p3 + p4) +
  plot_annotation(
    title   = "Energy News Sentiment: Three-Approach Comparison",
    caption = paste0("Synthetic data (", nrow(energy_news),
                     " articles). Replace with real RSS cache for ",
                     "research-grade results."),
    theme = theme(
      plot.title   = element_text(family = "serif", face = "bold",
                                   size = 14, colour = "#1e3a2f"),
      plot.caption = element_text(size = 8, colour = "#7a7a72")
    )
  )

Four-panel diagnostic for the energy sentiment comparison. Top left: weekly sentiment time series by approach. Top right: WTI price with sentiment overlay. Bottom left: score distributions. Bottom right: weekly article volume.

What this example teaches

The energy-price comparison illustrates several principles from the earlier sections in a concrete, applied setting:

Multi-lexicon comparison is essential. The three approaches produce scores on different scales (integer vs. proportion vs. continuous), with different distributional shapes, and with different sensitivity to negation. The correlation between them is often moderate, not high, which means your conclusions depend on your tool choice.

Topic segmentation changes the story. Aggregate sentiment across all energy articles masks the fact that supply-risk news is structurally negative while demand-recovery news is mixed. Segmenting by topic before aggregating produces more interpretable and potentially more predictive signals.

External validation is the gold standard. When you have an external criterion (price data), you can evaluate which approach performs best empirically, rather than arguing about which lexicon is theoretically superior. This is the same logic as the hand-coding validation in Section 8, but with a continuous external variable instead of human labels.

Modest correlations are normal. Sentiment-price correlations in energy markets are typically |r| = 0.05 to 0.20. If you find correlations above 0.30, investigate whether there is data leakage or whether the sentiment score is mechanically correlated with the price variable through shared keywords.

Challenges specific to environmental news

Environmental news has several properties that make sentiment analysis harder than in other domains:

Neutral-tone reporting. Wire services (AP, Reuters) follow style guides that enforce emotional flatness. This pushes lexicon scores toward zero and inflates the neutral class. Models trained on opinion text (tweets, reviews) systematically over-predict polarity on news wire text.

Attributed quotes. “Critics called the policy a disaster” reports someone else’s negative sentiment without the article endorsing it. Lexicon methods cannot distinguish the reporter’s voice from a quoted source. Quote detection (using dependency parsing via spacyr) can help.

Domain vocabulary. General-purpose lexicons mishandle domain-specific terms. “Stranded asset” is not a negative phrase about being abandoned; it is a precise financial term. “Loss and damage” is a policy concept, not an expression of grief. Domain-adapted models like ClimateBERT exist precisely to handle this vocabulary.

Negativity bias in coverage. Environmental news is structurally negative. Disasters, extinctions, and pollution events are news; stable ecosystems are not. Raw aggregate sentiment is therefore meaningless without a topic-conditional or publication-conditional baseline.

Multiple lexicons as a robustness check

Always compute sentiment under at least two independent lexicons and compare. If LSD2015 and AFINN produce time series that move in opposite directions on the same corpus, that tells you the signal is driven by lexicon composition, not by genuine shifts in news tone.

# Rank correlation between two lexicon time series
cor.test(
  monthly_lsd$mean_sentiment,
  monthly_afinn$mean_sentiment,
  method = "spearman"
)

Topic-conditional sentiment

Running an unsupervised topic model (LDA or BERTopic) on your corpus, then computing sentiment within each topic, reveals structure that aggregate indices hide. An overall “neutral” trend might mask simultaneous increases in negative wildfire sentiment and positive renewable-energy sentiment.

Data sources for environmental news

What’s next?

• Time-Series Analysis – decomposition, ARIMA, and forecasting for the temporal indices you build here
• Multiple Regression & GLMs – regress sentiment indices against external predictors (carbon prices, policy events)
• Data Visualization – polished ggplot2 figures for publication-quality sentiment plots
• Tidyverse Data Wrangling – the pipe-based data manipulation that underlies the entire workflow