Skip to content

Stochastic Frontier Analysis

0packages depend
7.7Kdownloads / year
–test coverage
13/13flavors OK

Overview

About

Provides a user-friendly framework for estimating a wide variety of cross-sectional and panel stochastic frontier models. Suitable for a broad range of applications, the implementation offers extensive flexibility in specification and estimation techniques.

Changes in 1.2.0

A feature release. In brief:

  • **A Wald test for whether the endogeneity correction was needed**, endogeneity_test(), from Hou, Ramalho and Roseta-Palma (2025). On an ivsfm() fit it tests H0: rho = 0 -- that the noise is uncorrelated with the reduced-form errors, so a plain sfm() fit would have been consistent. The null is INTERIOR, since rho lives in the open unit ball, so this is an ordinary chi-square and not the chi-bar-square mixture inefficiency_test() needs. Size 0.046 at a nominal 5% over 1000 replications; power 1.000 at rho = 0.4 and 0.8.
  • **ivsfm() now reports standard errors for rho.** They were NA, on the ground that rho = t/sqrt(1 + t't) has a non-diagonal Jacobian and a delta-method value computed as though it were diagonal would be wrong. The first half was right and the conclusion was not: the Jacobian is J = (I - rho rho')/s, so Var(rho) = J Var(t) J' exactly. Checked against a numerical Jacobian, and validated end-to-end by the size of the test above -- a wrong Jacobian would show up there as a mis-sized test. The full covariance travels with the fit as $vcov_rho, which the joint test needs.
Show all changes in 1.2.0Hide the rest of the changes in 1.2.0
  • **ttsfm(model_name = "TTHN") no longer fabricates a log-density it cannot compute.** The likelihood needs D = Phi2(x1, 0; rho1) - Phi2(x2, 0; rho2), a difference of two bivariate normal CDFs of the same order of magnitude, and in parts of the parameter space it cancels completely: at sigma_v = 0.3, sigma_u = 1, sigma_w = 0.2 on a 400-observation draw one observation's D comes back as exactly 0 while both CDFs are O(0.1). pmax(D, .Machine$double.xmin) did not repair that -- it invented log(D) = -708.4 for the observation, and that single observation moved the summed objective by **715 log-units** in a surface whose real curvature is a few units per 0.1 step in log sigma. A cliff that size is what L-BFGS-B reports as ABNORMAL_TERMINATION_IN_LNSRCH. D is now rejected when it falls below the rounding error of its own subtraction, which leaves genuine small probabilities alone: a D of 6e-94 where both CDFs are also ~1e-94 has full significance and is kept. Three smaller repairs in the same branch: the exponentials for sigma_v/sigma_u/sigma_w are clipped at EXP_CLIP_UPPER as the TTNE branch already clipped its own; rho is held strictly inside (-1, 1), which sigma_v small enough does reach; and the substitution of -sqrt(.Machine$double.xmax / n) -- about -9.5e152 at n = 200, 140 orders of magnitude past the finite penalty it sat beside -- is gone, replaced by that penalty. It had also mapped a +Inf log-density to a large NEGATIVE one, so an unusable draw could be rewarded. Measured on eight matched fits: one that previously failed outright now converges, and one that previously returned sigma_w = 0.044 against a truth of 1 now returns 1.047 **with a better likelihood**. The other six are unchanged. sigma_v still collapses toward zero on some draws; that is a separate problem and is not fixed.
  • **ttsfm() no longer throws away a fit because of the optimizer's exit code.** It was the only entry point calling stop() on ANY non-zero convergence code, including L-BFGS-B's code 52, ABNORMAL_TERMINATION_IN_LNSRCH, which it returns whenever its line search meets a discontinuity and usually *having improved* the objective. The guard was written when nothing checked opt$convergence at all; opt.optim() now rebuilds at the stage-2 point whenever the value or the Hessian is non-finite, so what was left of it only destroyed usable fits. Two of eight TTHN test fits returned code 52 with sound estimates where the stop() had aborted them.
  • **What skewness and kurtosis a pair of distributions can actually produce**, moment_range(), from Papadopoulos and Parmeter (2021). spec_test() already shipped their formal statistic; this is the range check that comes before it. Because the variance share R lies in (0,1), the composed error can never be more skewed than u itself -- Normal-Half-Normal is confined to gamma1 in (-0.995, 0), Normal-Exponential to (-2, 0) -- so residual moments outside that set refute the pair with no test statistic involved. That is their Type II failure, alongside the Type I failure of wrong skew. Five noise and five inefficiency distributions, including the two whose shape is a free parameter: the truncated normal sweeps gamma1(u) in (0, 2) and gamma2(u) in (-0.243, 6) and is the only u here whose composed error can be platykurtic under normal noise, while the gamma is unbounded in both and so can **never** suffer a Type II failure. The kurtosis range is not the interval between the two components: with leptokurtic noise it dips below both, to 0.674 for Laplace noise with a half-normal u against the half-normal's own 0.869. **It is a diagnostic, not a test, and the print method says so.** PP2021's Tables 1-5 replicate -- their population rows exactly, and Tables 2, 3 and 5 cell for cell -- and they show the range check over-rejects: with correctly specified data at n = 200 and SNR 3, the sample excess kurtosis leaves its theoretical range 26% of the time. Use spec_test() for the p-value.
  • **A frontier that tolerates the wrong skewness**, esfm(), from Hafner, Manner and Simar (2018). When the sample skewness of the residuals comes out positive -- a small-sample accident, not evidence against the model -- the classical MLE collapses: sigma_u is exactly zero, every firm is fully efficient and the fit says nothing. This model keeps the frontier and widens the distribution of u, with one parameter gamma whose SIZE is the scale of inefficiency and whose SIGN chooses the direction of skewness. The classical model is nested at gamma > 0. On a 50-observation wrongly skewed draw, sfm(model_name = "NHN") returns sigma_u = 0.0003 and mean efficiency 0.99976; esfm() returns gamma = -0.57 with the slopes still near truth. Because gamma = 0 is an **interior** point rather than a boundary, the LR test of "no inefficiency" is an ordinary chi2(1) -- not the chi-bar-square mixture inefficiency_test() needs for the classical model. Size over 2000 replications on the paper's design, at a nominal 5%: 0.178 / 0.082 / 0.054 / 0.061 / 0.058 for n = 50 / 100 / 200 / 500 / 1000, against the paper's 0.152 / 0.076 / 0.061 / 0.054 / 0.053. Badly oversized below n = 100, usable from n = 200 -- their own conclusion. **It is not a free improvement everywhere.** The paper reports smaller bias in mean efficiency than the classical model even when the population skewness has the correct sign, and that reproduces when inefficiency is not small: at gamma = 0.5 the bias is -0.0001 at n = 100 against the classical +0.021. But at gamma = 0.3, where 22% of samples are wrongly skewed, this model overshoots the other way (-0.040 against +0.036) and is the worse of the two at n = 50 (-0.080 against +0.044). A small true gamma lets it fit sizeable inefficiency of either sign to noise. Use it when the wrong skewness is the problem.
  • **Kim and Schmidt (2008) two-step test**, uhet_test(). The procedure applied work actually runs -- fit the frontier, regress the JLMS predictor on firm characteristics, test the slopes -- is invalid as usually performed, because u-hat is a *generated* dependent variable and the first-step estimation error enters the second-step variance. The correction is not distribution-free, so "NHN" and "NE" are implemented, the two cases the paper works through. What it buys: with z independent of x the corrected and naive tests agree, as the theory requires. With z strongly correlated with x (rho = 0.9) the naive standard error is **52% too large** at n = 800 and the naive test's size collapses to **0.000** at a nominal 5%, while the corrected standard error matches the true sampling spread to within half a percent (ratio 0.995) and keeps size 0.032. **Caveat, stated plainly:** the corrected test is conservative in these simulations. On the paper's own base case it rejects 0.024/0.036/0.028 at N = 200/500/1000 against a nominal 0.05, where the paper reports 0.0495 at N = 200. The *naive* branch does reproduce the paper's figure for the uncorrected test (0.023-0.029 against their 0.0237). The implementation was checked against equations (5)-(13), the first-step scores agree with an analytic half-normal score to 1e-11, gamma-hat is unbiased and Gaussian, and an oracle test using the true sampling standard deviation has correct size (0.045) -- so the gap is in the variance estimate, not in gamma-hat or the reference distribution. Two candidate explanations were tested and **both refuted**: the homoskedastic second-step variance (the paper's "OLS" against its "OLS-H") changes nothing, and neither does using the Hessian information rather than the outer product when forming r_i. The gap is unexplained rather than merely unverified. Read a rejection as meaningful and a non-rejection as weak evidence. Fits on the wrong-skew boundary are refused rather than answered: every E[u|eps] is then identical, so there is no estimated inefficiency for z to explain, and the information matrix is singular there (Waldman 1982).
  • **Chen and Wang (2012) centered-residuals moment test**, cw_test(). Unlike gof_test(), it never evaluates the composed-error density or distribution function: it compares the empirical characteristic function of the *centred* residuals with its theoretical value, which follows from the characteristic functions of the two components separately -- in closed form for many pairs whose convolution is not, the normal-gamma model being the standard example. Centring is what makes it work. A frontier's intercept is not identified separately from E[u], so the least-squares intercept estimates alpha - E[u]; centring cancels it, which is why the test is valid off an ordinary regression rather than requiring maximum likelihood. The variance corrects for both nuisances -- estimating the scales, and the centring itself. Omitting the second understates it and oversizes the test. Size at a nominal 5% over 500 replications, half-normal data: | test | n | tau=1 | tau=1.5 | tau=2 | |---|---|---|---|---| | cosine | 500 | 0.070 | 0.054 | 0.050 | | cosine | 2000 | 0.046 | 0.044 | 0.050 | | sine | 500 | 0.092 | 0.064 | 0.042 | | sine | 2000 | 0.068 | 0.062 | 0.064 | Power against exponential inefficiency at tau = 1 is 0.231 (n = 500) and 0.952 (n = 2000) for the cosine test against 0.108 and 0.706 for the sine test. Both of the authors' conclusions reproduce: the cosine test is the less sensitive to tau and the more powerful, so it is the default. **tau defaults to a single frequency**: passing several at once gives rejection rates of 0.09-0.15 at a nominal 5%, because the cosine moments at nearby frequencies are near-collinear and the covariance matrix is then near-singular. The authors report the same for their combined tau. A wrongly signed third moment leaves the moment equations with no admissible solution, and that is an error rather than a number.
  • **The Gauss-Legendre rule is now cached.** It depends on nothing but the node count and was rebuilt by an eigendecomposition on every call -- hundreds of times per cw_test(), since numDeriv differences a quadrature-based moment function. Caching also keeps eigen() out of forked workers: macOS's threaded BLAS is not fork-safe and mclapply() segfaulted inside it before this was added.
  • **Goodness-of-fit tests for the assumed inefficiency distribution**, gof_test(), from Wang, Amsler and Schmidt (2011). The distribution of u is the assumption in this model least often defended and least often tested, and it is testable: hold the normality of v fixed and the assumed u implies a distribution for the composed error, so rejecting that is rejecting the assumed u. Kolmogorov-Smirnov and Pearson chi-square, both on the composed error rather than on u-hat -- which is the same test, since u-hat = E[u|eps] is monotonic in eps, but far easier, and it avoids the trap the paper warns about: comparing the spread of u-hat with the assumed density of u is a mistake, not a diagnostic, because those are different distributions. Parameter estimation is not ignored. The default parametric bootstrap copies the estimation step exactly -- each replication draws from the fitted model, rebuilds the response, **refits**, forms its own residuals and recomputes the statistic. The "asymptotic" option gives a chi-square p-value that the authors note is conservative at the MLE, and **returns NA rather than a KS p-value**, because with parameters estimated the Kolmogorov distribution does not apply and Bai's (2003) transformation is not implemented. Size and power over 200 replications (half-normal data for size, exponential data fitted as half-normal for power), at nominal 0.10 / 0.05: | | n | KS | chi2 | |---|---|---|---| | size | 200 | 0.065 / 0.025 | 0.085 / 0.065 | | size | 800 | 0.070 / 0.040 | 0.065 / 0.010 | | power | 200 | 0.505 / 0.395 | 0.235 / 0.145 | | power | 800 | 0.985 / 0.965 | 0.810 / 0.655 | | power | 2000 | 1.000 / 1.000 | 0.990 / 0.985 | Both hold their size, mildly conservatively, and **KS dominates chi-square at every sample size** -- the paper's own conclusion, and the reason it is the better default.
  • **Coelli (1995)'s hypothesis tests, and his COLS standard errors.** New inefficiency_test() tests H0: no technical inefficiency, reporting the one-sided LR test, the naive LR test, the Wald ratio and the third-moment test side by side. The point of the paper is that **the two tests most often reported are the two with the wrong size**: H0 puts gamma on the boundary, so the LR statistic follows the Gourieroux-Holly-Monfort mixture 0.5 chi2(0) + 0.5 chi2(1) (5% critical value 2.71, not 3.84) and the Wald ratio is not asymptotically normal at all. Measured here over 400 replications per cell at a nominal 5%, with sigma_u = 0: | n | LR one-sided | LR naive | Wald | M3T | |---|---|---|---|---| | 100 | 0.060 | 0.030 | 0.230 | 0.055 | | 400 | 0.055 | 0.033 | 0.217 | 0.045 | | 800 | 0.060 | 0.035 | 0.182 | 0.062 | The one-sided LR and third-moment tests hold their size; the naive LR test rejects about half as often as it should; the Wald test rejects three to four times too often and does not improve with n. Under the alternative the one-sided LR test has the better power of the two valid tests at every sample size, which is Coelli's recommendation. NTN fits get the two-restriction mixture 0.25/0.5/0.25 over chi2(0,1,2).
  • **sfm(estimator = "cols") now reports Coelli's analytic standard errors.** It previously returned NA for sigv and sigu unless a bootstrap was requested, because the OLS errors are not appropriate for them: those parameters are non-linear functions of the residual moments m2 and m3, not regression coefficients. Coelli's Appendix 1 (A15-A16) gives the correct delta-method variance and it is now filled in for NHN; the intercept additionally carries the error in its E[u] shift. NE and NG keep NA and the bootstrap, because the moment inversion -- and therefore the variance -- is distribution-specific. Checked against the nonparametric bootstrap that was already there, on the same data: the analytic/bootstrap ratio for sigu is 1.047 at n = 200, 1.011 at n = 800 and 0.984 at n = 3200, which is an asymptotic approximation behaving as one should. Coelli's A13 inversion turned out to be the moment inversion the package already carried as .gtre_two_step() (equal to ten decimal places, not merely similar), and his A15 variance is the delta method whose GTRE_SEQ2 implementation is corrected elsewhere in this release -- **including the factor of 2 on the covariance term, which A15 independently confirms**, and A16's `mu6 - mu3^2 - 6 mu2 mu4 + 9 mu2^3

The NEWS section is longer than shown here.

From NEWS.md in the release

Install

Health

CRAN checks
13OK
Slowest check: 12.1 min · r-oldrel-windows-x86_64
Status changes
  • OK13 OK2026-08-25
  • WARNING12 OK, 1 WARNING2026-08-20
  • OK13 OK2026-08-04
  • NOTE12 OK, 1 NOTE2026-08-01
  • OK13 OK2026-06-09
Show 3 earlier changesHide 3 earlier changes
  • ERROR12 OK, 1 ERROR2026-06-08
  • WARNING12 OK, 1 WARNING2026-06-07
  • OK14 OK2026-03-10
Advisories

none known

Documentation

In release 1.2.0
READMEYes · 2,632 wordsVignettesYes · 1 · dynamicNEWSYes · 100% structured · mentions 1.2.0Citation fileYes · cites a publicationWebsiteYes
Help pages
47includes a package overview page
Examples that run
43%20 run in a plain check · 25 only when donttest examples are enabled (R CMD check --as-cran) · 1 never run
Exported functions with examples
100%
Documented parameters
84%
Return-value docs
100%
Help pages citing references
68%
Vignettes (1)

Downloads

7.7K
CRAN downloads in the past year
Rank #6,453 · ~21/day · ~640/mo
Daily download trend is not available in this view yet.
46830 days
1.3K90 days
7.7K1 year
Compare downloads with other packages →
Also on298 r2u2 autocran132 c2d4u

Repository

2Stars
1Forks
1Open issues
1Open PRs
0Releases
243Commits
1Commit identities
sfastochastic-frontierstochastic-frontier-analysis
License file reads as GPL-2.0 · 243 commits · Last activity 2026-10-16 · 0% stars, 30d · Default branch at 1.2.1, ahead of CRAN 1.2.0 (checked 2026-10-18)
Responsiveness
How these are sampled
PRs merged (all time)
96%
Median days to close, 50 most recently updated closed issues
0
Median days to close, 50 most recently updated closed PRs
0

Stars over time

2026-08-24 · 12026-09-16 · 2
PracticesBeta

3 development-tooling and community-health practices detected across 2 families in the upstream repository

Checks run against github.com/davidhbernstein/sfa on 2026-10-18.

Continuous integration (2)
GitHub ActionsR CMD check workflow

The R CMD check workflow runs on Linux, macOS and Windows, and on R-devel

CRAN release process (1)
cran-comments.md
How this is detected·Detection ruleset v3 (2026-09-26)
AI coding toolsBeta

This repository shows signs of AI coding tools.

  • Claude Codefirst seen 2026-08-20
    Pull request it wrote: opened by a project member
    Commits: at least 211 crediting it in the message
Models named in its commits
  • Claude Opus 5101 commits · Aug to Sep 2026
  • Claude Opus 5.56 commits · Sep 2026

Some of these counts are a minimum, because the search stopped before the oldest commits.

Counts are commits, not lines of code. Commits by a tool's accounts count once across all of its accounts. Commits crediting it are the largest count among its credit lines, a minimum.How this is detected·Detection ruleset v4 (updated 2026-09-27)

Dependencies

Declared dependencies
21 external dependencies (excludes base and recommended)

Installing sfa adds 27 packages (16 with compiled code, 14 MB of source tarballs). Heaviest: Rcpp 1.8 MB, collapse 1.4 MB, randtoolbox 1.4 MB.

Depends (1)
R >= 4.0.0
LinkingTo (0)
none
Enhances (0)
none
Reverse dependencies
0direct
0indirect

Nothing depends on this yet.

Datasets

People & History

People (3)
Author, Maintainer since 1.0.4 (2026) · added in 1.0.4
Author · added in 1.0.4
Author · added in 1.0.4

Previously maintained by Ariane Straub (1.0-0 to 1.0-1).

Listed in earlier versions (2)
under the supervision of Torsten Hothorn
no longer listed · 1.0-0 to 1.0-1
Ariane Straub
no longer listed · 1.0-0 to 1.0-1
Releases (5)

Pick two releases to compare their code metrics.

  • 1.2.0Latest
    record
    2026-09-07 · current release · diff ↗ ·
    added 33 exports
    lcsfm, selsfm, ivsfm, copsfm, efficiency_ci, skewness_test, marginal_effects, efficiency, meanefficiency, simulation_se, TIC, vuong, lcsfm_homogeneity, spec_test, spec_test_all, sfma, influence_sfa, hscore, hscore_select, calibrate_c, density_weights, pcomposed, dcomposed, inefficiency_test, pcomposed_model, composed_cdf, gof_test, cw_test, uhet_test, esfm and 3 more
  • 1.1.5
    record
    2026-08-23 · diff ↗ ·
    added 4 exports
    ttsfm, npsfm, sfa_diagnostics, psfm_bootstrap
  • R
    R 4.6.0 released · 2026-04-24
  • unarchivedReturned to CRAN
    2026-01-21
  • 1.0.4
    record
    2026-01-21 · diff ↗ ·
    removed 5 exports
    sfa, dgp, eff, u.sfa, te.eff.sfa
    added 7 exports
    zsfm, psfm, sfm, data_gen_p, data_gen_cs, summary.sfareg, print.sfareg
  • R
    R 4.5.0 released · 2025-04-11
  • R
    R 4.4.0 released · 2024-04-24
  • R
    R 4.3.0 released · 2023-04-21
  • archivedRemoved from CRAN
    2022-06-15
    check problems were not corrected despite reminders
  • R
    R 4.2.0 released · 2022-04-22
  • R
    R 4.1.0 released · 2021-05-18
  • R
    R 4.0.0 released · 2020-04-24
  • R
    R 3.6.0 released · 2019-04-26
  • R
    R 3.5.0 released · 2018-04-23
  • R
    R 3.4.0 released · 2017-04-21
  • R
    R 3.3.0 released · 2016-05-03
Show 9 earlier eventsHide 9 earlier events
  • R
    R 3.2.0 released · 2015-04-16
  • R
    R 3.1.0 released · 2014-04-10
  • 1.0-1
    record
    2014-01-06 · diff ↗
  • R
    R 3.0.0 released · 2013-04-03
  • R
    R 2.15.0 released · 2012-03-30
  • R
    R 2.14.0 released · 2011-10-31
  • R
    R 2.13.0 released · 2011-04-13
  • 1.0-0
    record
    2011-01-24
  • R
    R 2.12.0 released · 2010-10-15

Package metadata

First published
2011-01-24
Total releases
5 / 15 yrs
License
GPL (>= 2) OSI
Minimum R
≥ 4.0.0
Bundled data
391 KB / 4 files
Download size
1.0 MB

Cite

Cite this package

Run in R for the authors' entry:

citation("sfa")

The package provides a citation file, so citation() returns the authors' own entry.

Cite the R Observatory

For a number measured here: a download total, a coverage figure, an archival date.

APA

Balamuta, J. J. (2026). R Observatory: Metrics for sfa version 1.2.0 [Data set]. HJJB, LLC. Data release v2026-10-22. https://doi.org/10.5281/zenodo.21843040

From data release v2026-10-22, which the citation names so these numbers can be found later. More on citing and the projects behind them.

Report a problem with this page →

Privacy choices

These apply to this browser and are stored on this device only. Nothing about your choice is sent to us.

Read the privacy policy