GB generation share trends

NG ESO data

1 About

Code copyright (c) 2024 the author

License:

• CC-BY

To cite:

• Anderson, B, (2024) GB generation share trends: NG ESO data. https://dataknut.github.io/gridCarbon/gbGenMixTrends.html

Code:

• https://github.com/dataknut/gridCarbon

2 Introduction

Inspired by https://twitter.com/DrSimEvans/status/1508409309775994892

3 Data

This analysis uses the UK NG ESO generation mix data.

The data contains average (mean) MW generation per half hour by source including inter-connectors and storage.

As far as we can work out this data does not include distributed (i.e. non-grid connected) generation such as small scale wind, solar, hydro, biomass etc which is connected to the LV network. This means the ESO data is likely to underestimate total generation and potentially underestimate the proportion of total generation that is renewable. It is possible that this could be fixed using embedded wind & solar generation data from the demand update file.

The data also includes a half-hourly carbon intensity value in g CO2/kWh sourced from https://www.carbonintensity.org.uk/.

rmdParams$olderThan <- 7

If the data we have previously downloaded is more than 7` days old, re-download.

gbGenMixUrl <- "https://data.nationalgrideso.com/backend/dataset/88313ae5-94e4-4ddc-a790-593554d8c6b9/resource/f93d1835-75bc-43e5-84ad-12472b180a98/download/df_fuel_ckan.csv"

orig_DT <- get_gbGenMix(url = gbGenMixUrl, 
                      dataPath = repoParams$ukGridDataLoc,
                      olderThan = rmdParams$olderThan)

## Most recent date in local version of data is 2024-11-30 (1 days ago)... loading

orig_DT[, dv_dateTime := lubridate::as_datetime(DATETIME)] # proper date time
message("Original data range from: ", min(orig_DT$dv_dateTime))

## Original data range from: 2009-01-01

message("...to: ", max(orig_DT$dv_dateTime))

## ...to: 2024-11-30 21:00:00

# add derived variables used later ----
orig_DT[, dv_year := lubridate::year(dv_dateTime)]
orig_DT[, dv_date := lubridate::as_date(dv_dateTime)]
orig_DT[, dv_month := lubridate::month(dv_dateTime)]
orig_DT[, dv_hour := lubridate::hour(dv_dateTime)]
orig_DT[, dv_hms := hms::as_hms(dv_dateTime)]
# half-hours are the start of the half hours (we think)

orig_DT <- gridCarbon::add_season(orig_DT, 
                                 dateVar = "dv_dateTime",
                                 h = "N") # north

#message("Remove incomplete years to avoid weird things in plots.")
# remove incomplete days (can cause weired effects)
#gbGenMix_dt <- orig_DT[dv_year < 2023]
gbGenMix_dt <- orig_DT[dv_date < max(dv_date)]

message("Filtered data range from: ", min(gbGenMix_dt$dv_dateTime))

## Filtered data range from: 2009-01-01

message("...to: ", max(gbGenMix_dt$dv_dateTime))

## ...to: 2024-11-29 23:30:00

gbGenMix_dt <- gridCarbon::add_peakPeriod(gbGenMix_dt, 
                                      dateTime = "dv_dateTime")

# check coding
# table(gbGenMix_dt$dv_hour, gbGenMix_dt$dv_peak, useNA = "always")

rmdParams$plotCap <- paste0("Data: NG ESO Generation Mix ",
                   min(gbGenMix_dt$dv_date), " - ",
                   max(gbGenMix_dt$dv_date),
                            "\nPlot: @dataknut \nCode: https://github.com/dataknut/gridCarbon")

Note: * the data covers more years than we need * the data may contain partial years - BEWARE incomplete years in plots using annual totals or means across all months.

4 Recreating DrSimEvans’ plot (data from 2012 onwards)

This looks like daily data.

• Solar + wind (% of gen) vs coal + gas (& of gen))

# add together the % and totals we want (half-hourly)
gbGenMix_dt[, dv_coal_gas_pc := COAL_perc + GAS_perc]
gbGenMix_dt[, dv_solar_wind_pc := SOLAR_perc + WIND_perc]
gbGenMix_dt[, dv_coal_gas := COAL + GAS]
gbGenMix_dt[, dv_solar_wind := SOLAR + WIND]

# keep the vars we want for clarity
temp <- gbGenMix_dt[, .(dv_dateTime, dv_coal_gas, dv_solar_wind,
                    dv_coal_gas_pc, dv_solar_wind_pc, GENERATION)]

temp[, dv_date := lubridate::date(dv_dateTime)]

# aggregate to daily data for plotting
plotDT <- temp[,
               .(mean_dv_solar_wind_pc = mean(dv_solar_wind_pc),
                 mean_dv_coal_gas_pc = mean(dv_coal_gas_pc),
                 total_dv_coal_gas = sum(dv_coal_gas),
                 total_dv_solar_wind = sum(dv_solar_wind),
                 total_GENERATION = sum(GENERATION),
                 nObs = .N), # to check for days with < 48 half hours
               keyby = .(dv_date)
               ]
plotDT[, dv_year := lubridate::year(dv_date)] # for plots
plotDT[, total_dv_coal_gas_pc := total_dv_coal_gas/total_GENERATION] # daily %
plotDT[, total_dv_solar_wind_pc := total_dv_solar_wind/total_GENERATION]

message("Check for days with less than 48 hours - this will be truncated data")

## Check for days with less than 48 hours - this will be truncated data

table(plotDT$nObs)

## 
##   48 
## 5812

Figure 4.1 shows the mean half-hourly % generation by each type per day. This is slightly convoluted - it is the mean of the sum of the 48 daily half-hourly XXX_perc values in the original data where XXX is the generation type. Unfold the code above for clarity.

The smoothed curves are estimated for each year. The lines terminate at the maximum value for the year. I’m still trying to decide if they tell us anything useful.

ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = mean_dv_solar_wind_pc, 
                            y = mean_dv_coal_gas_pc,
                            colour = as.factor(dv_year),
                            alpha = dv_year)) +
  geom_point() +
  geom_smooth() +
  scale_colour_viridis_d(name = "Year") +
  guides(alpha = "none") +
  labs(x = "Solar & wind (mean % of half-hourly generation per day)",
       y = "Coal & gas (mean % of half-hourly generation per day)",
       caption = rmdParams$plotCap)

## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 4.1: Mean half-hourly % generation by each type per day

# save it
ggplot2::ggsave(filename = "meanGBrenewablesVsfossilHalfHourPC.png", 
                path = rmdParams$plotPath,
                height = 5)

## Saving 7 x 5 in image
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 4.2 shows the percentage of daily generation by type. This is less convoluted as it is the sum of generation per day for the two categories (solar + wind vs gas + coal) as a % of total daily generation.

Again the smoothed curve is estimated for each year.

ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = 100 * total_dv_solar_wind_pc, 
                            y = 100 * total_dv_coal_gas_pc,
                            colour = as.factor(dv_year),
                            alpha = dv_year)) +
  geom_point() +
  geom_smooth() +
  scale_colour_viridis_d(name = "Year") +
  guides(alpha = "none") +
  labs(x = "Solar & wind (% of total daily generation)",
       y = "Coal & gas (% of total daily generation)",
       caption = rmdParams$plotCap)

## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

Figure 4.2: Percentage of daily generation by type

ggplot2::ggsave(filename = "dailyGBpcGenMix.png", 
                path = rmdParams$plotPath,
                height = 5)

## Saving 7 x 5 in image
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

4.1 Half-hourly versions of the plot

Just cos we can… helpfully split into ‘peak’ and ‘off peak’ periods.

Peak period definitions:

• Morning 07:00 - 09:00
• Daytime 09:00 - 16:00
• Evening 16:00 - 21:00
• Night - all other times

Again the smoothed curve is estimated for each year (and demand period).

ggplot2::ggplot(gbGenMix_dt[dv_year > 2011], aes(x = dv_solar_wind_pc, 
                            y = dv_coal_gas_pc,
                            alpha = dv_year,
                            colour = as.factor(dv_year))) +
  geom_point() +
  facet_wrap(. ~ dv_peakPeriod) +
  geom_smooth() +
  scale_colour_viridis_d(name = "Year") +
  guides(alpha = "none") +
  labs(x = "Solar & wind (% of half-hourly generation)",
       y = "Coal & gas (% of half-hourly generation)",
       caption = rmdParams$plotCap)

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 4.3: Percentage of half-hourly generation by type

ggplot2::ggsave(filename = "halfHourlyPCgenByPeakPeriod.png", 
                path = rmdParams$plotPath,
                height = 5)

## Saving 7 x 5 in image
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

5 MWh generation trends

Values are in MW per half hour (presumably the mean over the half-hour) so divide by 2 to give MWh.

ref: https://data.nationalgrideso.com/carbon-intensity1/historic-generation-mix/r/historic_gb_generation_mix

dt_long <- melt(gbGenMix_dt, id.vars = c("dv_dateTime", "dv_year", "dv_hms", "dv_peakPeriod", "dv_hour"))

## Warning in melt.data.table(gbGenMix_dt, id.vars = c("dv_dateTime", "dv_year", : 'measure.vars' [DATETIME, GAS,
## COAL, NUCLEAR, ...] are not all of the same type. By order of hierarchy, the molten data value column will be
## of type 'character'. All measure variables not of type 'character' will be coerced too. Check DETAILS in
## ?melt.data.table for more on coercion.

dt_long[, dv_date := lubridate::date(dv_dateTime)]

plotDT <- dt_long[variable == "GAS" | variable == "COAL" |
                    variable == "NUCLEAR" | variable == "WIND" |
                    variable == "HYDRO" | variable == "IMPORTS" |
                    variable == "HYDRO" | variable == "BIOMASS" | 
                    variable == "OTHER" | variable == "SOLAR" |
                    variable == "STORAGE", 
                  .(dailyMWh = sum(as.numeric(value)/2)), # MW -> MWh
                  keyby = .(dv_year, dv_date, variable)]

# annual sum for cross-check 
t <- plotDT[, .(sum_TWh = sum(dailyMWh/1000000)), keyby = .(Year = dv_year)]
t[, Year := as.factor(Year)]
gridCarbon::make_flexTable(t, caption = "Total TWh per year for cross-check")

``` ```{=html}

(\#tab:longTrendsDaily)Total TWh per year for cross-check

Year	sum_TWh
2009	329.2
2010	335.2
2011	322.3
2012	323.8
2013	322.4
2014	309.9
2015	304.7
2016	303.1
2017	300.7
2018	298.4
2019	291.5
2020	275.5
2021	284.8
2022	289.7
2023	274.2
2024	253.6

# double check

Note that these values are lower then the ones found in https://assets.publishing.service.gov.uk/government/uploads/system/uploads/attachment_data/file/1043323/Regional_Electricity_Generation_and_Supply_2016-2020.pdf although they are of the same order of magnitude (200 - 300 TWh per year).

One reason for this is that Northern Ireland is excluded from the NG ESO data used here.

Figure 5.1 shows individual trends.

p <- ggplot2::ggplot(plotDT, aes(x = dv_date, y = dailyMWh/1000, colour = variable)) +
  geom_line() +
  facet_grid(variable ~ .) +
  scale_color_viridis_d(name = "Source") +
  labs(x = "Date",
       y = "GWh",
       caption = rmdParams$plotCap)
p

Figure 5.1: Trends in daily GWh generation 2012-2022

ggplot2::ggsave(filename = "dailyGenBySource.png", 
                path = rmdParams$plotPath,
                height = 5)

## Saving 7 x 5 in image

Figure 5.2 stacks them - spot COVID19 lockdown 2020…

ggplot2::ggplot(plotDT, aes(x = dv_date, y = dailyMWh/1000, fill = variable)) +
  geom_col(position = "stack") +
  scale_fill_viridis_d(name = "Source") +
  labs(x = "Date",
       y = "GWh",
       caption = rmdParams$plotCap)

Figure 5.2: Trends in daily GWh generation 2012-2022 (stacked)

ggplot2::ggsave(filename = "genTrendStackBySource.png", 
                path = rmdParams$plotPath
                )

## Saving 7 x 5 in image

5.1 Renewable %

Defined how?

we think renewable is wind + solar, low carbon includes nuclear - does this look likely?

# we think renewable is wind + solar, low carbon includes nuclear
gbGenMix_dt[, Year := lubridate::year(dv_date)]
ggplot2::ggplot(gbGenMix_dt, aes(x = RENEWABLE_perc, 
                                   y = LOW_CARBON_perc,
                                   colour = NUCLEAR_perc)) +
  geom_point() +
  facet_wrap(Year ~ .)

make_NgesoPlots <- function(dt, 
                            var = "RENEWABLE_perc", # defaults
                            lowColour = "green", 
                            highColour = "red", 
                            scaleLab = "Change me!", 
                            minMax = "max"){
  res <- list() # for the results
  res$tile <- ggplot2::ggplot(dt, aes(x = dv_date, 
                                      y = dv_hms, fill = get(var))) +
    geom_tile() +
    theme(legend.position = "bottom") +
    scale_x_date(date_labels="%b %Y",date_breaks  ="12 month") +
    theme(axis.text.x = element_text(angle = 60, hjust = 1)) +
    scale_fill_continuous(name = scaleLab,
                          high = highColour,
                          low = lowColour) +
    labs(x = "Date",
         y = "Half-hour",
         caption = rmdParams$plotCap)
  
  # if(minMax = "max"){
  #   res$tile <- res$tile + 
  #     scale_fill_continuous(name = paste0("Maximum ", scaleLab)
  #                           )
  # }
  #   if(minMax = "min"){
  #   res$tile <- res$tile + 
  #     scale_fill_continuous(name = paste0("Half-hourly ", scaleLab)
  #                           )
  #   }
  
  # line plot
  plotDT <- dt[,.(max = max(get(var)),
                  mean = mean(get(var)),
                  min = min(get(var))
  ),
  keyby = .(dv_date)]
  
  
  
  p <- ggplot2::ggplot(plotDT, aes(x = dv_date)) +
    theme(legend.position = "bottom") +
    scale_x_date(date_labels="%b %Y",date_breaks  ="12 month") +
    theme(axis.text.x = element_text(angle = 60, hjust = 1))
  
  if(minMax == "max"){
    yh <- max(plotDT$max)
    p <- p + geom_line(aes(y = max, colour = max)) +
      scale_color_continuous(name = paste0("Maximum half-hourly ", scaleLab),
                             high = highColour,
                             low = lowColour) +
      geom_smooth(aes(y = max), colour = "grey") +
      geom_hline(yintercept = yh, colour = highColour) +
      annotate("text", x = mean(plotDT$dv_date), 
                      y = yh*1.05, label = paste0("Maximum: ", round(yh)))
    p <- p + labs(x = "Date",
                  y = paste0("Maximum half-hourly ", scaleLab),
                  caption = rmdParams$plotCap
    )
  }
  
  if(minMax == "min"){
    yh <- min(plotDT$min)
    p <- p + geom_line(aes(y = min, colour = min)) +
      scale_color_continuous(name = paste0("Minimum half-hourly ", scaleLab),
                             high = highColour,
                             low = lowColour) +
      geom_smooth(aes(y = min), colour = "grey") + 
  geom_hline(yintercept = yh, colour = lowColour) +
      annotate("text", x = mean(plotDT$dv_date), 
                      y = yh*1.05, label = paste0("Minimum: ", round(yh)))
    p <- p + labs(x = "Date",
                  y = paste0("Minimum half-hourly ", scaleLab),
                  caption = rmdParams$plotCap
    )
  }
  res$line <- p
  
  return(res)
}

rpc <- make_NgesoPlots(gbGenMix_dt, 
                       var = "RENEWABLE_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Renewable generation (%)",
                       minMax = "max" # max line or min line?
)
rpc$tile

Figure 5.3: Half-hourly % renewables

rpc$line

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 5.4: Half-hourly % renewables

rpc <- make_NgesoPlots(gbGenMix_dt[Year > 2021], 
                       var = "RENEWABLE_perc",
                       lowColour = "red",
                       highColour = "green",
                       scaleLab = "Renewable generation (%)",
                       minMax = "max"
)
rpc$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="3 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 5.5: Recent trends

rpc$line +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="3 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 5.6: Recent trends

5.2 Wind %

#plotDT <- gbGenMix_dt[, .()]
ci <- make_NgesoPlots(gbGenMix_dt, 
                      var = "WIND_perc",
                      lowColour = "red",
                      highColour = "green",
                      scaleLab = "Wind %",
                      minMax = "max"
)
ci$tile

Figure 5.7: Wind %

ci$line

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 5.8: Wind %

wpc <- gbGenMix_dt[, .(mean_pc = mean(WIND_perc)), keyby = .(dv_hms, Year)]

ggplot2::ggplot(wpc, aes(x = dv_hms, colour = Year, group = Year, y = mean_pc)) +
  geom_line() +
  labs(x = "Time of day",
       y = "% wind generation")

Figure 5.9: Wind %

ci <- make_NgesoPlots(gbGenMix_dt[Year > 2021], 
                      var = "WIND_perc",
                      lowColour = "red",
                      highColour = "green",
                      scaleLab = "Wind %",
                      minMax = "max"
)
ci$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="3 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 5.10: Recent wind %

ci$line +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="3 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 5.11: Recent wind %

6 Carbon intensity trends

Now let’s use the CARBON_INTENSITY variable.

6.1 Overall trends

Seasons:

• Summer: Jun - Aug
• Autumn: Sept - Nov
• Winter: Dec - Feb
• Spring: Mar - May

plotDT <- gbGenMix_dt[, .(mean_CI = mean(CARBON_INTENSITY)), keyby = .(dv_year, dv_date, dv_season)]

ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = dv_date, y = mean_CI, colour = dv_season)) +
  geom_point() +
  scale_color_viridis_d(name = "Season") +
  theme(legend.position = "bottom") +
  labs(x = "Date",
       y = "Mean g CO2e/kWh",
       caption = rmdParams$plotCap)

Figure 6.1: Mean half-hourly carbon intensity per day 2012-2022

ggplot2::ggsave(filename = "meanCiTrendBySeason.png", 
                path = rmdParams$plotPath
                )

## Saving 7 x 5 in image

Re-draw as a boxplot of mean daily CI by month - plotted at month start. 6.2 shows outliers nicely.

plotDT <- gbGenMix_dt[, .(mean_CI = mean(CARBON_INTENSITY)), 
                     keyby = .(dv_year, 
                               dv_date,
                               dv_season)]
plotDT[, dv_month := lubridate::floor_date(dv_date, unit = "months"),]
ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = dv_month, 
                                            y = mean_CI, 
                                            group = dv_month,
                                            colour = dv_season)) +
  geom_boxplot() +
  scale_color_viridis_d(name = "Season") +
  theme(legend.position = "bottom") +
  labs(x = "Month",
       y = "Mean g CO2e/kWh",
       caption = rmdParams$plotCap)

Figure 6.2: Mean daily CI boxplots per month

ggplot2::ggsave(filename = "meanCiTrendByMonthSeason.png", 
                path = rmdParams$plotPath
                )

## Saving 7 x 5 in image

Figure 6.3 shows half-hourly carbon intensity over time.

#plotDT <- gbGenMix_dt[, .()]
ci <- make_NgesoPlots(gbGenMix_dt, 
                      var = "CARBON_INTENSITY",
                      lowColour = "green",
                      highColour = "brown",
                      scaleLab = "Carbon intensity (g CO2/MW)",
                      minMax = "min"
)
ci$tile

Figure 6.3: Half-hourly carbon intensity over time

ci$line

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

ci <- make_NgesoPlots(gbGenMix_dt[Year > 2021], 
                      var = "CARBON_INTENSITY",
                      lowColour = "green",
                      highColour = "brown",
                      scaleLab = "Carbon intensity (g CO2/MW)",
                      minMax = "min"
)
ci$tile +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="3 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.

Figure 6.4: Recent

ci$line +
  scale_x_date(date_labels="%b-%Y",date_breaks  ="3 month")

## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

Figure 6.5: Recent

6.2 Annual by peak period

Mean carbon intensity per year and season within peak period.

Figure 6.6 suggests evening peak periods still have slightly higher carbon intensity and the shape of the reduction curves differ by season although rather less by period. Interestingly the sustained reduction in carbon intensity in Summer has leveled off.

plotDT <- gbGenMix_dt[, .(mean_CI = mean(CARBON_INTENSITY)), keyby = .(dv_year, dv_peakPeriod, dv_season)]

ggplot2::ggplot(plotDT, aes(x = dv_year, y = mean_CI, colour = dv_peakPeriod)) +
  geom_line() +
  scale_color_viridis_d(name = "Peak period") +
  facet_wrap(. ~ dv_season) +
  labs(x = "Year",
       y = "Mean g CO2e/kWh",
       caption = rmdParams$plotCap)

Figure 6.6: Mean half-hourly carbon intensity by peak period and season 2012-2022

ggplot2::ggsave(filename = "annualMeanCIByPeak.png", 
                path = rmdParams$plotPath
                )

## Saving 7 x 5 in image

7 Renewables and overall generation

Do we see a relationship between renewables generating and peak demand? This will be mediated by the way the electricity market works.

We may find wind curtailment (not visible here) at low demand periods where nuclear can’t be shut off.

First, what is the general average shape of carbon intensity and renewable generation?

Figure 7.1 shows that although the mean half-hourly carbon intensity had fallen over time (with the effect of solar in summer particularly noticeable), the morning and evening peaks are still relatively more carbon intense as demand overtakes the available renewable supply.

plotDT <- gbGenMix_dt[dv_year > 2011, .(mean_ci = mean(CARBON_INTENSITY),
                          mean_renewables_MW = mean(RENEWABLE),
                          mean_renewables_pc = mean(RENEWABLE_perc)),
                          keyby = .(dv_year, dv_hms, dv_season)]

ggplot2::ggplot(plotDT, aes(x = dv_hms, y = mean_ci, 
                            alpha = dv_year,
                            colour = dv_year,
                            group = dv_year)) +
  geom_line() +
  scale_alpha_continuous(name="Year") +
  scale_color_continuous(name = "Year",
                         low = "grey", 
                         high = "#3CBAC6") + # UoS from Marine palette
  theme(axis.text.x = element_text(angle = 45, hjust=1)) +
  facet_wrap(. ~ dv_season) +
  labs(x = "Time of day",
       y = "Mean half-hourly carbon intensity",
       caption = rmdParams$plotCap)

Figure 7.1: Mean half-hourly carbon intensity by year and season

ggplot2::ggplot(plotDT, aes(x = dv_hms, y = mean_renewables_MW/1000, 
                            colour = dv_year, 
                            alpha = dv_year,
                            group = dv_year)) +
  geom_line() +
  scale_alpha_continuous(name = "Year") +
  scale_color_continuous(name = "Year",
                         low = "grey", 
                         high = "#3CBAC6") + # UoS from Marine palette
  facet_grid(dv_season ~ .) +
    labs(x = "Time of day",
       y = "Mean renewables (GW)",
       caption = rmdParams$plotCap)

Figure 7.2: Mean half-hourly renewable generation by year and season

ggplot2::ggplot(plotDT, aes(x = dv_hms, y = mean_renewables_pc, 
                            colour = dv_year, 
                            alpha = dv_year,
                            group = dv_year)) +
  geom_line() +
  scale_alpha_continuous(name = "Year") +
  scale_color_continuous(name = "Year",
                         low = "grey", 
                         high = "#3CBAC6") + # UoS from Marine palette
  facet_grid(dv_season ~ .) +
    labs(x = "Time of day",
       y = "Mean % renewables (%)",
       caption = rmdParams$plotCap)

Figure 7.3: Mean half-hourly renewable generation by year and season

ggplot2::ggplot(gbGenMix_dt, aes(x = RENEWABLE, y = GENERATION)) +
  geom_point() +
  facet_wrap(. ~ dv_year) +
  theme(axis.text.x = element_text(angle = 45, hjust=1))

7.4 shows the rise of both wind and solar (mid-day peak)…

Beware incomplete yearss

plotDT <- gbGenMix_dt[, .(mean_renewables = mean(RENEWABLE),
                          mean_solar = mean(SOLAR),
                 mean_generation = mean(GENERATION)),
             keyby = .(dv_year, dv_hms, dv_peakPeriod, dv_season)]

ggplot2::ggplot(plotDT[dv_year > 2016], aes(x = dv_hms, 
                                         colour = dv_peakPeriod,
                                         alpha = dv_year,
                                         group = dv_year)) +
  geom_line(aes(y = mean_renewables/1000)) +
  scale_alpha_continuous(name = "Year") +
  scale_color_discrete(name = "Peak period") +
  theme(axis.text.x = element_text(angle = 45, hjust=1)) +
  facet_wrap(. ~ dv_season) +
  labs(x = "Time of day",
       y = "Mean renewables (GW)",
       caption = rmdParams$plotCap)

Figure 7.4: Trends in mean half-hourly renewable generation by time of day and season

What would have happened if we increased solar generation in 2022 by a factor of 10? We would have needed some storage in Spring & Summer…

ggplot2::ggplot(plotDT[dv_year == 2022], aes(x = dv_hms,
                                             colour = dv_peakPeriod,
                                             group = dv_peakPeriod)) +
  geom_point(aes(y = mean_generation/1000)) +
  geom_line(aes(y = mean_solar/1000 * 10), linetype = "dotdash") +
  scale_color_discrete(name = "Peak period") +
    theme(axis.text.x = element_text(angle = 45, hjust=1)) +
  facet_wrap(. ~ dv_season) +
  labs(x = "Time of day",
       y = "Mean generation (GW)",
       caption = rmdParams$plotCap)

Figure 7.5: Comparing mean total generation & 10 * solar generation by time of day for 2020

YMMV on 7.6

plotDT <- gbGenMix_dt[, .(mean_renewables = mean(RENEWABLE),
                 mean_generation = mean(GENERATION)),
             keyby = .(dv_year, dv_hms, dv_peakPeriod, dv_season)]

ggplot2::ggplot(plotDT[dv_year > 2016], aes(x = mean_generation/1000 , y = mean_renewables/1000,
                                         colour = dv_peakPeriod)) +
  geom_point() +
  scale_color_discrete(name = "Period") +
  facet_grid(dv_season ~ dv_year) +
  labs(x = "Mean total generation (GW)",
       y = "Mean renewables (GW)",
       caption = rmdParams$plotCap)

Figure 7.6: Mean renewables vs Mean total generation

8 Summary

That’s it.

You might want to look at recent academic research on this topic:

• (Staffell 2017)
• (Staffell and Pfenninger 2018)

9 Annex

9.1 Data descriptors

9.1.1 NGESO gen mix data

skimr::skim(gbGenMix_dt)

Table 9.1: Data summary

Name	gbGenMix_dt
Number of rows	278976
Number of columns	46
Key	NULL
_______________________
Column type frequency:
Date	1
difftime	2
factor	2
numeric	39
POSIXct	2
________________________
Group variables	None

Variable type: Date

skim_variable	n_missing	complete_rate	min	max	median	n_unique
dv_date	0	1	2009-01-01	2024-11-29	2016-12-15	5812

Variable type: difftime

skim_variable	n_missing	complete_rate	min	max	median	n_unique
dv_hms	0	1	0 secs	84600 secs	42300 secs	48
hms	0	1	0 secs	84600 secs	42300 secs	48

Variable type: factor

skim_variable	n_missing	complete_rate	ordered	n_unique	top_counts
dv_season	0	1	FALSE	4	Spr: 70656, Sum: 70656, Aut: 69840, Win: 67824
dv_peakPeriod	0	1	FALSE	5	Ear: 81368, Day: 81368, Eve: 46496, Lat: 46496

Variable type: numeric

skim_variable	n_missing	complete_rate	mean	sd	p0	p25	p50	p75	p100	hist
GAS	0	1	12237.95	5420.92	566.0	7841.00	12321.0	16461.00	27472.0	▃▇▇▅▁
COAL	0	1	5929.39	6580.38	0.0	294.00	2745.0	11207.00	26044.0	▇▂▂▁▁
NUCLEAR	0	1	6456.55	1426.89	2065.0	5281.00	6660.0	7657.00	9342.0	▁▃▆▇▅
WIND	0	1	4716.47	4418.14	1.0	1239.00	3299.0	6873.00	21998.0	▇▃▂▁▁
HYDRO	0	1	393.54	245.12	0.0	189.00	362.0	560.00	1403.0	▇▇▅▁▁
IMPORTS	0	1	2365.81	1526.91	0.0	1324.00	2354.0	3024.00	9148.0	▆▇▂▁▁
BIOMASS	0	1	847.15	1053.00	0.0	0.00	0.0	1776.00	3328.0	▇▁▂▂▁
OTHER	0	1	458.98	601.08	0.0	44.00	161.0	722.00	2462.0	▇▂▁▁▁
SOLAR	0	1	862.74	1744.22	0.0	0.00	0.0	792.00	11472.0	▇▁▁▁▁
STORAGE	0	1	279.99	360.41	0.0	0.00	158.0	420.00	2660.0	▇▁▁▁▁
GENERATION	0	1	34548.56	7225.54	18341.0	28906.75	34316.0	39598.00	59577.0	▃▇▇▂▁
CARBON_INTENSITY	0	1	303.01	146.77	19.0	183.00	275.0	440.00	644.0	▃▇▅▆▂
LOW_CARBON	0	1	13276.45	4990.50	4626.0	9326.00	12075.0	16385.00	36297.0	▇▇▃▁▁
ZERO_CARBON	0	1	12429.30	4486.63	3799.0	9010.00	11290.0	15122.00	35428.0	▇▇▃▁▁
RENEWABLE	0	1	5972.75	4999.34	2.0	1839.00	4619.0	8862.00	30049.0	▇▃▂▁▁
FOSSIL	0	1	18167.33	9154.82	566.0	11205.00	17224.5	24326.00	49096.0	▅▇▅▂▁
GAS_perc	0	1	34.89	13.00	1.8	25.00	35.8	45.10	72.7	▂▆▇▆▁
COAL_perc	0	1	15.72	16.34	0.0	1.00	7.9	31.30	60.6	▇▂▂▂▁
NUCLEAR_perc	0	1	19.33	5.40	5.0	15.40	18.7	22.50	42.3	▁▇▅▁▁
WIND_perc	0	1	14.42	13.83	0.0	3.60	9.8	21.10	69.5	▇▃▂▁▁
HYDRO_perc	0	1	1.12	0.65	0.0	0.60	1.1	1.60	4.4	▇▇▃▁▁
IMPORTS_perc	0	1	7.32	5.13	0.0	3.80	6.9	9.90	36.6	▇▆▁▁▁
BIOMASS_perc	0	1	2.65	3.33	0.0	0.00	0.0	5.60	16.0	▇▃▂▁▁
OTHER_perc	0	1	1.37	1.82	0.0	0.10	0.5	2.00	10.3	▇▂▁▁▁
SOLAR_perc	0	1	2.45	5.05	0.0	0.00	0.0	2.20	38.3	▇▁▁▁▁
STORAGE_perc	0	1	0.72	0.90	0.0	0.00	0.5	1.10	7.9	▇▁▁▁▁
GENERATION_perc	0	1	100.00	0.00	100.0	100.00	100.0	100.00	100.0	▁▁▇▁▁
LOW_CARBON_perc	0	1	39.98	16.36	10.7	26.90	37.1	50.50	91.6	▅▇▅▂▁
ZERO_CARBON_perc	0	1	37.32	14.70	10.7	25.90	34.5	46.10	88.9	▅▇▅▂▁
RENEWABLE_perc	0	1	17.99	15.21	0.0	5.40	13.8	26.70	74.6	▇▅▂▁▁
FOSSIL_perc	0	1	50.62	18.88	2.4	37.30	51.9	65.80	88.0	▂▅▇▇▅
dv_year	0	1	2016.46	4.59	2009.0	2012.00	2016.0	2020.00	2024.0	▇▆▆▆▆
dv_month	0	1	6.49	3.43	1.0	4.00	7.0	9.00	12.0	▇▆▆▆▇
dv_hour	0	1	11.50	6.92	0.0	5.75	11.5	17.25	23.0	▇▇▆▇▇
dv_coal_gas_pc	0	1	50.62	18.88	2.4	37.30	51.9	65.80	88.0	▂▅▇▇▅
dv_solar_wind_pc	0	1	16.87	15.10	0.0	4.40	12.6	25.50	74.3	▇▃▂▁▁
dv_coal_gas	0	1	18167.33	9154.82	566.0	11205.00	17224.5	24326.00	49096.0	▅▇▅▂▁
dv_solar_wind	0	1	5579.21	4957.54	1.0	1482.00	4190.0	8427.00	29545.0	▇▃▂▁▁
Year	0	1	2016.46	4.59	2009.0	2012.00	2016.0	2020.00	2024.0	▇▆▆▆▆

Variable type: POSIXct

skim_variable	n_missing	complete_rate	min	max	median	n_unique
DATETIME	0	1	2009-01-01	2024-11-29 23:30:00	2016-12-15 23:45:00	278976
dv_dateTime	0	1	2009-01-01	2024-11-29 23:30:00	2016-12-15 23:45:00	278976

9.1.2 Check definitions

It is not clear from https://data.nationalgrideso.com/carbon-intensity1/historic-generation-mix/r/historic_gb_generation_mix how the following are defined:

• LOW_CARBON
• ZERO_CARBON
• RENEWABLE
• FOSSIL

test2021 <- gbGenMix_dt[dv_year == 2021]

#test2021[, ba_LOW_CARBON := ]

9.2 CI model

Model the half-hourly carbon intensity - just for fun but so we can use it elsewhere (if we assume the same relationships hold!).

We use just coal, gas, wind, hydro & solar as they also appear in the NZ gen mix that we’d like to apply to the model to.

9.2.1 MWh generation models

Model relationship between gas, coal, wind, hydro and solar generation and carbon intensity as a linear model. Would also be interesting to test a neural network model…

mod1 <- lm(CARBON_INTENSITY ~ GAS + COAL + WIND + HYDRO + SOLAR, data = gbGenMix_dt)

ggstats::ggcoef_model(mod1, exponentiate = FALSE)

The results are pretty much what we’d expect…

As above but add year as a co-variate

mod2 <- lm(CARBON_INTENSITY ~ GAS + COAL + WIND + HYDRO + SOLAR + dv_year, data = gbGenMix_dt)

summary(mod2)

## 
## Call:
## lm(formula = CARBON_INTENSITY ~ GAS + COAL + WIND + HYDRO + SOLAR + 
##     dv_year, data = gbGenMix_dt)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -136.56  -16.71    1.57   17.31  117.28 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  8.418e+03  4.471e+01  188.29   <2e-16 ***
## GAS          2.584e-03  1.167e-05  221.40   <2e-16 ***
## COAL         1.506e-02  1.364e-05 1103.82   <2e-16 ***
## WIND        -7.868e-03  1.875e-05 -419.58   <2e-16 ***
## HYDRO       -9.198e-03  2.535e-04  -36.29   <2e-16 ***
## SOLAR       -6.505e-03  3.413e-05 -190.62   <2e-16 ***
## dv_year     -4.061e+00  2.217e-02 -183.19   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 29.04 on 278969 degrees of freedom
## Multiple R-squared:  0.9609, Adjusted R-squared:  0.9609 
## F-statistic: 1.141e+06 on 6 and 278969 DF,  p-value: < 2.2e-16

ggstats::ggcoef_model(mod2)

9.2.2 % generation model

Try a % model - a % generation model is ‘normalised’ so with a lot of ceteris paribus we can apply it to other generation datasets assuming all the same relationships and carbon intensities are true. Which is really a bit fanciful…

9.2.2.1 All years (% model)

modpc <- lm(CARBON_INTENSITY ~ GAS_perc + COAL_perc + WIND_perc + HYDRO_perc + SOLAR_perc, data = gbGenMix_dt)

ggstats::ggcoef_model(modpc)

Interestingly hydro has a +ve % coeff although not when modelled as MWh. Hydro must be correlating with higher carbon generation to meet peaks?

9.2.2.2 2023 (% model)

Just 2023…

modpc2 <- lm(CARBON_INTENSITY ~ GAS_perc + COAL_perc + WIND_perc + HYDRO_perc + SOLAR_perc, data = gbGenMix_dt[dv_year == 2023])

ggstats::ggcoef_model(modpc2)

summary(modpc2)

## 
## Call:
## lm(formula = CARBON_INTENSITY ~ GAS_perc + COAL_perc + WIND_perc + 
##     HYDRO_perc + SOLAR_perc, data = gbGenMix_dt[dv_year == 2023])
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -17.0379  -4.4880  -0.3434   4.0446  26.6888 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 41.144997   0.433840   94.84   <2e-16 ***
## GAS_perc     3.573011   0.006787  526.48   <2e-16 ***
## COAL_perc   10.046020   0.041835  240.13   <2e-16 ***
## WIND_perc   -0.450460   0.006333  -71.12   <2e-16 ***
## HYDRO_perc   1.151356   0.070509   16.33   <2e-16 ***
## SOLAR_perc  -0.333893   0.009031  -36.97   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.076 on 17514 degrees of freedom
## Multiple R-squared:  0.9907, Adjusted R-squared:  0.9907 
## F-statistic: 3.713e+05 on 5 and 17514 DF,  p-value: < 2.2e-16

Table 9.2 reports full estimates so we can re-use them.

modpc2_df <- broom::tidy(modpc2, conf.int = TRUE)
knitr::kable(modpc2_df,
             digits = 3, 
             caption = "% generation model for half hourly carbon intentity (2023 only)")

Table 9.2: % generation model for half hourly carbon intentity (2023 only)

term	estimate	std.error	statistic	p.value	conf.low	conf.high
(Intercept)	41.145	0.434	94.839	0	40.295	41.995
GAS_perc	3.573	0.007	526.480	0	3.560	3.586
COAL_perc	10.046	0.042	240.134	0	9.964	10.128
WIND_perc	-0.450	0.006	-71.125	0	-0.463	-0.438
HYDRO_perc	1.151	0.071	16.329	0	1.013	1.290
SOLAR_perc	-0.334	0.009	-36.972	0	-0.352	-0.316

ggplot2::ggplot(modpc2_df, aes(x = term, y = estimate)) + 
  geom_point() + 
  geom_errorbar(aes(ymin=conf.low, 
                    ymax=conf.high),
                width = 0.2) +
  labs(title = "Model results",
       x = 'Variable',
       y = 'Coefficient',
       caption = paste0("Error bars = 95% CI")) +
  coord_flip() # rotate for legibility

Note that in both models hydro has a positive effect on CI which is unexpected.

9.3 R environment

Packages etc:

• base R (R Core Team 2016)
• bookdown (Xie 2016a)
• data.table (Dowle et al. 2015)
• ggplot2 (Wickham 2009)
• here (Müller 2017)
• hms (Müller 2018)
• knitr (Xie 2016b)
• lubridate (Grolemund and Wickham 2011)
• rmarkdown (Allaire et al. 2018)

References

Allaire, JJ, Yihui Xie, Jonathan McPherson, Javier Luraschi, Kevin Ushey, Aron Atkins, Hadley Wickham, Joe Cheng, and Winston Chang. 2018. Rmarkdown: Dynamic Documents for r. https://CRAN.R-project.org/package=rmarkdown.

Dowle, M, A Srinivasan, T Short, S Lianoglou with contributions from R Saporta, and E Antonyan. 2015. Data.table: Extension of Data.frame. https://CRAN.R-project.org/package=data.table.

Grolemund, Garrett, and Hadley Wickham. 2011. “Dates and Times Made Easy with lubridate.” Journal of Statistical Software 40 (3): 1–25. http://www.jstatsoft.org/v40/i03/.

Müller, Kirill. 2017. Here: A Simpler Way to Find Your Files. https://CRAN.R-project.org/package=here.

———. 2018. Hms: Pretty Time of Day. https://CRAN.R-project.org/package=hms.

R Core Team. 2016. R: A Language and Environment for Statistical Computing. Vienna, Austria: R Foundation for Statistical Computing. https://www.R-project.org/.

Staffell, Iain. 2017. “Measuring the Progress and Impacts of Decarbonising British Electricity.” Energy Policy 102 (March): 463–75. https://doi.org/10.1016/j.enpol.2016.12.037.

Staffell, Iain, and Stefan Pfenninger. 2018. “The Increasing Impact of Weather on Electricity Supply and Demand.” Energy 145 (February): 65–78. https://doi.org/10.1016/j.energy.2017.12.051.

Wickham, Hadley. 2009. Ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag New York. http://ggplot2.org.

Xie, Yihui. 2016a. Bookdown: Authoring Books and Technical Documents with R Markdown. Boca Raton, Florida: Chapman; Hall/CRC. https://github.com/rstudio/bookdown.

———. 2016b. Knitr: A General-Purpose Package for Dynamic Report Generation in r. https://CRAN.R-project.org/package=knitr.