Introduction to pgt

library(pgt)

The materials-balance principle

A pollution-generating technology turns material inputs into a good output and an unavoidable bad output. The pollutant content of the inputs has to end up somewhere: embodied in the product, removed by abatement, or emitted. pgt builds efficiency models on that accounting identity, \[u'x_l - v y_l \ge b_l,\] where \(x_l\) are the material inputs of unit \(l\), \(u\) their pollutant flow coefficients (for example tonnes of CO2 potential per unit of input), \(y_l\) the good output, \(v\) the pollutant retained in the product, and \(b_l\) the emissions. The gap \(u'x_l - v y_l - b_l\) is the pollutant that leaves the process other than through the measured emission. A note on terminology: the mass balance is an identity when every exit is measured; the inequality form above allows for unmeasured abatement, and when the abatement output \(a_l\) is observed the account closes as the equality \(u'x_l - v y_l = b_l + a_l\). Following Rødseth (2025), the documentation calls both forms the materials-balance identity.

The principle itself is old: Ayres and Kneese (1969) built the materials-balance view of production and externalities, Lauwers (2009) made the case for carrying it into frontier-based eco-efficiency models, and Dakpo, Jeanneaux and Latruffe (2016) survey how pollution-generating technologies have been modelled in nonparametric benchmarking since.

Existing R packages handle undesirable outputs by transforming the data (the Seiford and Zhu 2002 translation in deaR) or by imposing weak disposability and directional distances (nonparaeff, Benchmarking). Neither approach carries a materials-balance account. pgt audits the account before estimation, lets it enter the estimators as a cap, an equality, a price system or an objective valuation depending on the model, and puts the competing axiom systems behind one interface.

Building a technology

The constructor bundles the data with the flow coefficients. The shipped steeldemo data is a synthetic panel of steel plants on two production routes, with emissions in CO2 units, an observed abatement output (CO2 captured by carbon capture and storage, capture for utilisation, or slag mineralisation at the plants that run them) and the CO2 potential of the energy the capture consumes as a pollution-control input.

data(steeldemo)
tech <- pgt_tech(
  x = steeldemo[, c("coal_coke", "other_fuel", "raw_material", "flux",
                    "capture_energy")],
  y = steeldemo$production,
  b = steeldemo$emissions,
  a = steeldemo$captured,
  v = 0.01467,
  x_abate = "capture_energy",
  group = steeldemo$route,
  id = steeldemo$plant
)
tech
#> Pollution-generating technology
#>   DMUs: 180;  inputs: 5;  good outputs: 1;  bad outputs: 1
#>   inputs: coal_coke, other_fuel, raw_material, flux, capture_energy
#>   pollution-control inputs: capture_energy
#>   material flow coefficients: u = [1, 1, 1, 1, 1], v = [0.01467]
#>   abatement output: observed
#>   groups: Integrated (105), Minimill (75)
#>   materials balance: satisfied for all DMUs

The inputs already sit in pollutant-potential units, so the default \(u = 1\) applies to each; v = 0.01467 is the carbon retained per unit of crude steel (the synthetic generator’s assumption; see ?steeldemo for the derivation). captured fills the abatement slot a, and x_abate marks capture_energy as the pollution-control input, which the five-component decomposition below frees in its own stage while every estimator treats it like any other input.

Auditing the identity before estimation

mb_check() reports the closure gap per DMU and flags the accounts where more pollutant leaves than enters. The audit is the data-quality check behind the materials-balance interpretation of every model, and it gates the models that constrain the projection by the DMU’s own account: the input-fixed benchmark "wgd_input_fixed" and the directional model "fdmo".

mb <- mb_check(tech)
head(as.data.frame(mb))
#>     id      group potential retained       b gap rel_gap violated a closure
#> 1 P001 Integrated   5350250 33605.04 5316645   0       0    FALSE 0       0
#> 2 P001 Integrated   4528353 29414.52 4498938   0       0    FALSE 0       0
#> 3 P001 Integrated   4234863 28338.29 4206524   0       0    FALSE 0       0
#> 4 P002 Integrated   3629719 28435.22 3601283   0       0    FALSE 0       0
#> 5 P002 Integrated   3896048 30063.50 3865984   0       0    FALSE 0       0
#> 6 P002 Integrated   3933817 30179.46 3903637   0       0    FALSE 0       0
attr(mb, "n_violations")
#> [1] 0

With the abatement output observed, the audit also reports the closure of the full identity \(u'x - v y = b + a\) (the closure column), which is exact on this panel: every tonne of CO2 potential is accounted for by the product, the emission or the captured stream.

An infeasible "wgd_input_fixed" programme always belongs to a DMU with a strictly negative gap, so the audit locates any estimation problem in advance; such DMUs return NA scores and pgt() warns. The weak-G-disposability programme of Eq. 6 is always self-feasible, so its scores lie in \((0, 1]\) on any data.

Estimating environmental efficiency

The weak-G-disposability model of Rødseth (2025, Eq. 6) minimises emissions over peer mixes that hold the good output at the DMU’s level, pricing any extra good output at the DMU’s retained-content coefficient \(v\); the inputs are decision variables (the programmes are stated in full in vignette("models", "pgt")). The score is \(b^*_l / b_l\).

fit <- pgt(tech, model = "wgd")
summary(fit)
#> pgt fit: model = wgd, returns = vrs, peers = all
#>   pollutant: b
#>   DMUs: 180   failed LPs: 0
#> 
#> Efficiency (b*/b):
#>     0%    25%    50%    75%   100% 
#> 0.1056 0.1449 0.4402 0.6720 1.0000 
#> 
#> By group:
#>       group   n n_na median   mean median_dual_output
#>  Integrated 105    0 0.1566 0.2543             0.6921
#>    Minimill  75    0 0.6489 0.6607             0.2417

The envelope model drops the retained-content term (the \(v = 0\) case) and reduces the programme to the convex lower envelope of the \((y, b)\) scatter.

fit_env <- pgt(tech, model = "envelope")
summary(fit_env)
#> pgt fit: model = envelope, returns = vrs, peers = all
#>   pollutant: b
#>   DMUs: 180   failed LPs: 0
#> 
#> Efficiency (b*/b):
#>     0%    25%    50%    75%   100% 
#> 0.1056 0.1449 0.4402 0.6720 1.0000 
#> 
#> By group:
#>       group   n n_na median   mean median_dual_output
#>  Integrated 105    0 0.1566 0.2543             0.6921
#>    Minimill  75    0 0.6489 0.6606             0.2417

The directional representation of Rødseth (2025, Eq. 14) needs the closed account and holds abatement at the DMU’s own level; it expands the good output and contracts the bad output jointly, so its scores are in output units rather than ratios, and the materials-balance row ties the two together as \(\theta_b = v\,\theta_y\).

fit_dir <- pgt(tech, model = "fdmo")
summary(fit_dir)
#> pgt fit: model = fdmo, returns = vrs, peers = all
#>   pollutant: b
#>   DMUs: 180   failed LPs: 0
#> 
#> Gross inefficiency (theta_y + theta_b):
#>        0%       25%       50%       75%      100% 
#>       0.0       0.0  396565.8  706019.7 1554459.8 
#> 
#> By group:
#>       group   n n_na median   mean median_good median_bad
#>  Integrated 105    0 607786 602699      598999       8787
#>    Minimill  75    0 135782 221199      133818       1963

Shadow prices and abatement cost

The output-constraint dual \(\eta_l = \partial b^*/\partial y\) (returned as dual_output) measures the marginal emission content of the good output along the frontier. mac_curve() turns it into a marginal abatement cost curve: reducing the bad output by one unit costs \(1 / \eta_l\) units of the good output, or \(p / \eta_l\) at output price \(p\). Read the curve with care about which margin is which: each DMU’s plotted abatement is its distance to the frontier, \(b_l - b^*_l\), which the model itself prices at zero output loss, while the mac value is the marginal cost of abating beyond the DMU’s frontier point. The curve therefore ranks DMUs by the shadow price at their projection; the area under it is not a total-cost estimate.

head(shadow_prices(fit))
#>     id      group       b   b_star dual_output
#> 1 P001 Integrated 5316645 602998.0   0.2774621
#> 2 P001 Integrated 4498938 523740.4   0.2774621
#> 3 P001 Integrated 4206524 503385.0   0.2774621
#> 4 P002 Integrated 3601283 505218.3   0.2774621
#> 5 P002 Integrated 3865984 536014.9   0.2774621
#> 6 P002 Integrated 3903637 538208.1   0.2774621
mac <- mac_curve(fit, price = 550)
plot(mac)

Step curve of marginal abatement cost against cumulative abatement potential

Decomposing across production routes

With a technology group, pgt_decompose() splits environmental efficiency into a within-group and a technology-gap component. The envelope decomposition is an exact identity, Total \(=\) WR \(\times\) TGR, where WR is within-route reallocation and TGR the technology-gap ratio of the metafrontier tradition (Battese, Rao and O’Donnell 2004; O’Donnell, Rao and Battese 2008).

dec <- pgt_decompose(tech, type = "envelope")
summary(dec)
#> pgt decomposition: type = envelope, returns = vrs, DMUs = 180
#>   components: WR x TGR (product = total)
#> 
#> Group medians (components, total, attribution shares):
#>       group   n n_na    WR   TGR total share_WR share_TGR
#>  Integrated 105    0 0.183 0.990 0.157    0.997   0.00304
#>    Minimill  75    0 0.653 0.971 0.649    0.960   0.03963
#> 
#> (share_* are medians of DMU-level shares; they need not sum to 1.
#>  n_na counts DMUs with any NA component.)

The five-component source decomposition of Rødseth (2025, Eq. 11) traces each score to production technical efficiency, coefficient quality, production allocative efficiency, abatement technical efficiency and abatement allocative efficiency; it needs no groups, and the stage programmes are stated in vignette("models", "pgt"). The abatement components have something to measure here because the panel observes captured and marks capture_energy as the control input.

dec5 <- pgt_decompose(tech, type = "rodseth")
summary(dec5)
#> pgt decomposition: type = rodseth, returns = vrs, DMUs = 180
#>   components: te_production x quality x ae_production x te_abatement x ae_abatement (product = total)
#> 
#> Group medians (components, total, attribution shares):
#>       group   n n_na te_production quality ae_production te_abatement
#>  Integrated 105    0             1       1         0.425            1
#>    Minimill  75    0             1       1         0.861            1
#>  ae_abatement total share_te_production share_quality share_ae_production
#>         0.756 0.157           -1.46e-13             0               0.842
#>         0.932 0.649            2.08e-13             0               0.765
#>  share_te_abatement share_ae_abatement
#>            2.09e-17             0.0507
#>            0.00e+00             0.1381
#> 
#> (share_* are medians of DMU-level shares; they need not sum to 1.
#>  n_na counts DMUs with any NA component.)
comps <- c("te_production", "quality", "ae_production",
           "te_abatement", "ae_abatement", "total")
med <- aggregate(dec5$results[comps],
                 list(abatement = steeldemo$abatement_tech), median)
med[comps] <- round(med[comps], 3)
med
#>        abatement te_production quality ae_production te_abatement ae_abatement
#> 1           none         1.000       1         0.542        1.000        0.785
#> 2            CCS         1.000       1         1.000        0.911        0.990
#> 3            CCU         0.964       1         0.349        0.986        0.955
#> 4 mineralisation         1.000       1         0.992        0.566        0.850
#>   total
#> 1 0.287
#> 2 0.693
#> 3 0.145
#> 4 0.450

Two features of the table follow from the stage programmes. Quality is 1 throughout because the coefficients are homogeneous. Production technical efficiency is 1 at the plants without capture and below 1 at many capturing plants: the first stage holds the captured tonnage at the plant’s own level only as a floor, so a mix of peers that captures more from the same inputs reaches a lower controlled emission. Production allocative efficiency, the stage that frees the material inputs, carries most of the gap. The abatement components separate the plants by what they capture: abatement technical efficiency falls below 1 only where capture is active, where it measures how far the plant’s captured tonnage falls short of the best capture observed at no more capture energy, and abatement allocative efficiency, the stage that frees the capture-energy input, is lowest at the plants that capture nothing, where it measures the emission reduction available from adopting capture on the peers’ terms.

A worked example on real data

The frontier package ships a balanced panel of Philippine rice farms with a physical fertiliser input (NPK, kilograms of active ingredient) and a physical output (PROD, tonnes of paddy). Treating the nitrogen carried by the fertiliser as the pollutant potential gives a nitrogen materials balance. The nitrogen fraction below is illustrative (urea-equivalent, about 0.46 of the active ingredient).

data("riceProdPhil", package = "frontier")
d8 <- riceProdPhil[riceProdPhil$YEARDUM == 8, ]
uN <- 0.46
rice <- pgt_tech(
  x = as.matrix(d8[, c("AREA", "LABOR", "NPK", "OTHER")]),
  y = d8$PROD,
  b = uN * d8$NPK,
  u = c(AREA = 0, LABOR = 0, NPK = uN, OTHER = 0),
  v = 0,
  id = as.character(d8$FMERCODE)
)
rice_fit <- pgt(rice, model = "wgd")
summary(rice_fit)
#> pgt fit: model = wgd, returns = vrs, peers = all
#>   pollutant: b
#>   DMUs: 43   failed LPs: 0
#> 
#> Efficiency (b*/b):
#>     0%    25%    50%    75%   100% 
#> 0.1341 0.3777 0.5148 0.7235 1.0000

Here environmental efficiency is the scope to cut fertiliser nitrogen while holding rice output at its observed level. One caveat on this construction: because the bad output is defined as the full nitrogen potential of a single input (\(b = u_N \cdot\) NPK), every farm’s account closes with equality by construction and mb_check() passes trivially; the example illustrates the workflow, not the audit.

The other vignettes state the estimating programmes in full (vignette("models", "pgt")), reproduce the published-result replications (vignette("replication", "pgt")), compare the axiom systems (vignette("comparing-axioms", "pgt")), demonstrate multi-pollutant technologies and heterogeneous coefficients (vignette("multiple-pollutants", "pgt")) and cover productivity measurement with inference (vignette("productivity", "pgt")).

References