Skip to content

Add Conditional Value at Risk (CVAR) metric #100

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
abdelrahman-ayad wants to merge 8 commits into
base: main
Choose a base branch
from
Open

Add Conditional Value at Risk (CVAR) metric #100

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
bd3c07a
initial CVAR functions commit
abdelrahman-ayad Jan 30, 2026
b1f105a
added normalized CVAR (NCVAR)
abdelrahman-ayad Feb 2, 2026
e633074
formulate normalized CVAR based on CVAR struct
abdelrahman-ayad Feb 2, 2026
fac883c
broadcast CVAR and NCVAR - to check
abdelrahman-ayad Feb 2, 2026
37ab22c
broadcast CVAR and NCVAR - to check
abdelrahman-ayad Feb 2, 2026
2091f49
add tests for CVAR and NCVAR
abdelrahman-ayad Feb 2, 2026
04e6f25
CVAR and NCVAR shortfall and metrics test
abdelrahman-ayad Feb 2, 2026
3547a87
broadcast CVAR and NCVAR functions
abdelrahman-ayad Feb 2, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
43 changes: 34 additions & 9 deletions PRAS.jl/test/runtests.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,15 +1,40 @@
using PRAS
using Test

sys = PRAS.rts_gmlc()
@testset "ShortfallResult" begin
sys = PRAS.rts_gmlc()

sf, = assess(sys, SequentialMonteCarlo(samples=100), Shortfall())
sf, = assess(sys, SequentialMonteCarlo(samples=100), Shortfall())

eue = EUE(sf)
lole = LOLE(sf)
neue = NEUE(sf)
eue = EUE(sf)
lole = LOLE(sf)
neue = NEUE(sf)

@test val(eue) isa Float64
@test stderror(eue) isa Float64
@test val(neue) isa Float64
@test stderror(neue) isa Float64
@test val(eue) isa Float64
@test stderror(eue) isa Float64
@test val(neue) isa Float64
@test stderror(neue) isa Float64
end

@testset "ShortfallSamplesResult" begin
sys = PRAS.rts_gmlc()

sf, = assess(sys, SequentialMonteCarlo(samples=100), ShortfallSamples())

eue = EUE(sf)
lole = LOLE(sf)
neue = NEUE(sf)

alpha = 0.95
cvar = CVAR(sf, alpha)
ncvar = NCVAR(sf, cvar)

@test val(eue) isa Float64
@test stderror(eue) isa Float64
@test val(neue) isa Float64
@test stderror(neue) isa Float64
@test val(cvar) isa Float64
@test stderror(cvar) isa Float64
@test val(ncvar) isa Float64
@test stderror(ncvar) isa Float64
end
18 changes: 16 additions & 2 deletions PRASCore.jl/src/Results/Results.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@ import Base: broadcastable, getindex, merge!
import OnlineStats: Series
import OnlineStatsBase: EqualWeight, Mean, Variance, value
import Printf: @sprintf
import StatsBase: mean, std, stderror
import StatsBase: mean, std, stderror, quantile

import ..Systems: SystemModel, ZonedDateTime, Period,
PowerUnit, EnergyUnit, conversionfactor,
Expand All @@ -13,7 +13,7 @@ export

# Metrics
ReliabilityMetric, LOLE, EUE, NEUE,
val, stderror,
val, stderror, CVAR, NCVAR,

# Result specifications
Shortfall, ShortfallSamples,
Expand Down Expand Up @@ -80,6 +80,20 @@ NEUE(x::AbstractShortfallResult, r::AbstractString, ::Colon) =
NEUE(x::AbstractShortfallResult, ::Colon, ::Colon) =
NEUE.(x, x.regions.names, permutedims(x.timestamps))

CVAR(x::AbstractShortfallResult, alpha::Float64, ::Colon, t::ZonedDateTime) =
CVAR.(x, alpha, x.regions.names, t)

CVAR(x::AbstractShortfallResult, alpha::Float64, r::AbstractString, ::Colon) =
CVAR.(x, alpha, r, x.timestamps)

CVAR(x::AbstractShortfallResult, alpha::Float64, ::Colon, ::Colon) =
CVAR.(x, alpha, x.regions.names, permutedims(x.timestamps))

NCVAR(x::AbstractShortfallResult, alpha::Float64, r::AbstractString, ::Colon) =
NCVAR.(x, alpha, r, x.timestamps)

NCVAR(x::AbstractShortfallResult, alpha::Float64, ::Colon, ::Colon) =
NCVAR.(x, alpha, x.regions.names, permutedims(x.timestamps))
include("Shortfall.jl")
include("ShortfallSamples.jl")

Expand Down
83 changes: 83 additions & 0 deletions PRASCore.jl/src/Results/ShortfallSamples.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -185,6 +185,89 @@ function NEUE(x::ShortfallSamplesResult, r::AbstractString)

end

function CVAR(x::ShortfallSamplesResult{N,L,T,P,E}, alpha::Float64) where {N,L,T,P,E}
estimate = x[]
tail_losses = estimate[estimate .>= quantile(estimate, alpha)]

cvar = if !isempty(tail_losses)
MeanEstimate(tail_losses)
else
MeanEstimate(0.)
end

return CVAR{N,L,T,E}(cvar, alpha)

end

function CVAR(x::ShortfallSamplesResult{N,L,T,P,E}, alpha::Float64, r::AbstractString) where {N,L,T,P,E}
estimate = x[r]
tail_losses = estimate[estimate .>= quantile(estimate, alpha)]

cvar = if !isempty(tail_losses)
MeanEstimate(tail_losses)
else
MeanEstimate(0.)
end

return CVAR{N,L,T,E}(cvar, alpha)

end

function CVAR(x::ShortfallSamplesResult{N,L,T,P,E}, alpha::Float64, t::ZonedDateTime) where {N,L,T,P,E}
estimate = x[t]
tail_losses = estimate[estimate .>= quantile(estimate, alpha)]

cvar = if !isempty(tail_losses)
MeanEstimate(tail_losses)
else
MeanEstimate(0.)
end

return CVAR{N,L,T,E}(cvar, alpha)

end

function CVAR(x::ShortfallSamplesResult{N,L,T,P,E}, alpha::Float64, r::AbstractString, t::ZonedDateTime) where {N,L,T,P,E}
estimate = x[r, t]
tail_losses = estimate[estimate .>= quantile(estimate, alpha)]

cvar = if !isempty(tail_losses)
MeanEstimate(tail_losses)
else
MeanEstimate(0.)
end

return CVAR{N,L,T,E}(cvar, alpha)

end

function NCVAR(x::ShortfallSamplesResult{N,L,T,P}, cvar::CVAR) where {N,L,T,P}
demand = sum(x.regions.load)

ncvar = if demand > 0
div(cvar.cvar, demand/1e6)
else
MeanEstimate(0.)
end

return NCVAR(ncvar, cvar.alpha)

end

function NCVAR(x::ShortfallSamplesResult{N,L,T,P}, cvar::CVAR, r::AbstractString) where {N,L,T,P}
i_r = findfirstunique(x.regions.names, r)
demand = sum(x.regions.load[i_r, :])

ncvar = if demand > 0
div(cvar.cvar, demand/1e6)
else
MeanEstimate(0.)
end

return NCVAR(ncvar, cvar.alpha)

end

function finalize(
acc::ShortfallSamplesAccumulator{S},
system::SystemModel{N,L,T,P,E},
Expand Down
62 changes: 62 additions & 0 deletions PRASCore.jl/src/Results/metrics.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -164,3 +164,65 @@ function Base.show(io::IO, x::NEUE)
print(io, "NEUE = ", x.neue, " ppm")

end

"""
CVAR

`CVAR` reports conditional value at risk of shortfalls.

Contains both the estimated value itself as well as the standard error
of that estimate, which can be extracted with `val` and `stderror`,
respectively.
"""
struct CVAR{N,L,T<:Period,E<:EnergyUnit} <: ReliabilityMetric

cvar::MeanEstimate
alpha::Float64

function CVAR{N,L,T,E}(cvar::MeanEstimate, alpha::Float64) where {N,L,T<:Period,E<:EnergyUnit}
val(cvar) >= 0 || throw(DomainError(
"$val is not a valid CVAR"))
new{N,L,T,E}(cvar, alpha)
end

end

val(x::CVAR) = val(x.cvar)
stderror(x::CVAR) = stderror(x.cvar)

function Base.show(io::IO, x::CVAR{N,L,T,E}) where {N,L,T,E}

print(io, "CVAR@$(x.alpha) = ", x.cvar, " ",
unitsymbol(E), "/", N*L == 1 ? "" : N*L, unitsymbol(T))

end

"""
NCVAR

`NCVAR` reports normalized conditional value at risk of shortfalls.

Contains both the estimated value itself as well as the standard error
of that estimate, which can be extracted with `val` and `stderror`,
respectively.
"""
struct NCVAR <: ReliabilityMetric

ncvar::MeanEstimate
alpha::Float64

function NCVAR(ncvar::MeanEstimate, alpha::Float64)
val(ncvar) >= 0 || throw(DomainError(
"$val is not a valid NCVAR"))
new(ncvar, alpha)
end

end

val(x::NCVAR) = val(x.ncvar)
stderror(x::NCVAR) = stderror(x.ncvar)

function Base.show(io::IO, x::NCVAR)
print(io, "NCVAR@$(x.alpha) = ", x.ncvar, " ppm")

end
24 changes: 24 additions & 0 deletions PRASCore.jl/test/Results/metrics.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -72,4 +72,28 @@

end

@testset "CVAR" begin

cvar1 = CVAR{2,1,Hour,MWh}(MeanEstimate(1.2), 0.95)
@test string(cvar1) == "CVAR@0.95 = 1.20000 MWh/2h"

cvar2 = CVAR{1,2,Year,GWh}(MeanEstimate(17.2, 1.3), 0.95)
@test string(cvar2) == "CVAR@0.95 = 17±1 GWh/2y"

@test_throws DomainError CVAR{1,1,Hour,MWh}(MeanEstimate(-1.2), 0.95)

end

@testset "NCVAR" begin

ncvar1 = NCVAR(MeanEstimate(1.2), 0.95)
@test string(ncvar1) == "NCVAR@0.95 = 1.20000 ppm"

ncvar2 = NCVAR(MeanEstimate(17.2, 1.3), 0.95)
@test string(ncvar2) == "NCVAR@0.95 = 17±1 ppm"

@test_throws DomainError NCVAR(MeanEstimate(-1.2), 0.95)

end

end
40 changes: 40 additions & 0 deletions PRASCore.jl/test/Results/shortfall.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -100,6 +100,7 @@ end
N = DD.nperiods
r, r_idx, r_bad = DD.testresource, DD.testresource_idx, DD.notaresource
t, t_idx, t_bad = DD.testperiod, DD.testperiod_idx, DD.notaperiod
alpha = 0.95

result = PRASCore.Results.ShortfallSamplesResult{N,1,Hour,MW,MWh,ShortfallSamples}(
Regions{N,MW}(DD.resourcenames, DD.resource_vals), DD.periods, DD.d)
Expand All @@ -123,6 +124,16 @@ end
@test val(neue) ≈ mean(result[]) / load*1e6
@test stderror(neue) ≈ std(result[]) / sqrt(DD.nsamples) / load*1e6

cvar = CVAR(result, alpha)
estimate = result[];
tail_losses = estimate[estimate .>= quantile(estimate, alpha)];
@test val(cvar) ≈ mean(tail_losses)
@test stderror(cvar) ≈ std(tail_losses) / sqrt(length(tail_losses))

ncvar = NCVAR(result, cvar)
@test val(ncvar) ≈ val(cvar) / load*1e6
@test stderror(ncvar) ≈ stderror(cvar) / load*1e6

# Region-specific

@test length(result[r]) == DD.nsamples
Expand All @@ -142,10 +153,22 @@ end
@test val(region_neue) ≈ mean(result[r]) / load*1e6
@test stderror(region_neue) ≈ std(result[r]) / sqrt(DD.nsamples) / load*1e6

region_cvar = CVAR(result, alpha, r)
region_estimate = result[r];
region_tail_losses = region_estimate[region_estimate .>= quantile(region_estimate, alpha)];
@test val(region_cvar) ≈ mean(region_tail_losses)
@test stderror(region_cvar) ≈ std(region_tail_losses) / sqrt(length(region_tail_losses))

region_ncvar = NCVAR(result, region_cvar, r)
@test val(region_ncvar) ≈ val(region_cvar) / load*1e6
@test stderror(region_ncvar) ≈ stderror(region_cvar) / load*1e6

@test_throws BoundsError result[r_bad]
@test_throws BoundsError LOLE(result, r_bad)
@test_throws BoundsError EUE(result, r_bad)
@test_throws BoundsError NEUE(result, r_bad)
@test_throws BoundsError CVAR(result, alpha, r_bad)
@test_throws BoundsError NCVAR(result, region_cvar, r_bad)

# Period-specific

Expand All @@ -161,9 +184,16 @@ end
@test val(period_eue) ≈ mean(result[t])
@test stderror(period_eue) ≈ std(result[t]) / sqrt(DD.nsamples)

period_cvar = CVAR(result, alpha, t)
period_estimate = result[t];
period_tail_losses = period_estimate[period_estimate .>= quantile(period_estimate, alpha)];
@test val(period_cvar) ≈ mean(period_tail_losses)
@test stderror(period_cvar) ≈ std(period_tail_losses) / sqrt(length(period_tail_losses))

@test_throws BoundsError result[t_bad]
@test_throws BoundsError LOLE(result, t_bad)
@test_throws BoundsError EUE(result, t_bad)
@test_throws BoundsError CVAR(result, alpha, t_bad)

# Region + period-specific

Expand All @@ -180,6 +210,12 @@ end
@test val(regionperiod_eue) ≈ mean(result[r, t])
@test stderror(regionperiod_eue) ≈ std(result[r, t]) / sqrt(DD.nsamples)

regionperiod_cvar = CVAR(result, alpha, r, t)
regionperiod_estimate = result[r, t];
regionperiod_tail_losses = regionperiod_estimate[regionperiod_estimate .>= quantile(regionperiod_estimate, alpha)];
@test val(regionperiod_cvar) ≈ mean(regionperiod_tail_losses)
@test stderror(regionperiod_cvar) ≈ std(regionperiod_tail_losses) / sqrt(length(regionperiod_tail_losses))

@test_throws BoundsError result[r, t_bad]
@test_throws BoundsError result[r_bad, t]
@test_throws BoundsError result[r_bad, t_bad]
Expand All @@ -192,4 +228,8 @@ end
@test_throws BoundsError EUE(result, r_bad, t)
@test_throws BoundsError EUE(result, r_bad, t_bad)

@test_throws BoundsError CVAR(result, alpha, r, t_bad)
@test_throws BoundsError CVAR(result, alpha, r_bad, t)
@test_throws BoundsError CVAR(result, alpha, r_bad, t_bad)

end
Loading