Reference
Exported
ArviZ.rcParams
— ConstantrcParams
Dictionary to contain ArviZ default parameters, with validation when setting items.
Examples
julia> rcParams["plot.backend"]
"matplotlib"
julia> rcParams["plot.backend"] = "bokeh"
"bokeh"
julia> rcParams["plot.backend"]
"bokeh"
ArviZ.InferenceData
— TypeInferenceData(::PyObject)
InferenceData(; kwargs...)
Loose wrapper around arviz.InferenceData
, which is a container for inference data storage using xarray.
InferenceData
can be constructed either from an arviz.InferenceData
or from multiple Dataset
s assigned to groups specified as kwargs
.
Instead of directly creating an InferenceData
, use the exported from_xyz
functions or convert_to_inference_data
.
ArviZ.autocorr
— FunctionSee documentation for arviz.autocorr
.
ArviZ.autocov
— FunctionSee documentation for arviz.autocov
.
ArviZ.bfmi
— FunctionSee documentation for arviz.bfmi
.
ArviZ.compare
— FunctionSee documentation for arviz.compare
.
ArviZ.concat
— FunctionSee documentation for arviz.concat
.
ArviZ.concat!
— Functionconcat!(data1::InferenceData, data::InferenceData...; kwargs...) -> InferenceData
In-place version of concat
, where data1
is modified to contain the concatenation of data
and args
. See concat
for a description of kwargs
.
ArviZ.convert_to_inference_data
— FunctionSee documentation for arviz.convert_to_inference_data
.
ArviZ.convert_to_inference_data
— Methodconvert_to_inference_data(obj::Chains; group = :posterior, kwargs...) -> InferenceData
Convert the chains obj
to an InferenceData
with the specified group
.
Remaining kwargs
are forwarded to from_mcmcchains
.
ArviZ.convert_to_inference_data
— Methodconvert_to_inference_data(obj::NamedTuple; kwargs...) -> InferenceData
convert_to_inference_data(obj::Vector{<:NamedTuple}; kwargs...) -> InferenceData
convert_to_inference_data(obj::Matrix{<:NamedTuple}; kwargs...) -> InferenceData
convert_to_inference_data(obj::Vector{Vector{<:NamedTuple}}; kwargs...) -> InferenceData
Convert obj
to an InferenceData
. See from_namedtuple
for a description of obj
possibilities and kwargs
.
ArviZ.convert_to_inference_data
— Methodconvert_to_inference_data(
obj::SampleChains.AbstractChain;
group=:posterior,
kwargs...,
) -> InferenceData
convert_to_inference_data(
obj::SampleChains.AbstractChain;
group=:posterior,
kwargs...,
) -> InferenceData
Convert the chains obj
to an InferenceData
with the specified group
.
Remaining kwargs
are forwarded to from_samplechains
.
ArviZ.ess
— FunctionSee documentation for arviz.ess
.
ArviZ.extract_dataset
— FunctionSee documentation for arviz.extract_dataset
.
ArviZ.from_cmdstan
— FunctionSee documentation for arviz.from_cmdstan
.
ArviZ.from_cmdstan
— Methodfrom_cmdstan(posterior::Chains; kwargs...) -> InferenceData
Call from_mcmcchains
on output of CmdStan
.
ArviZ.from_dict
— FunctionSee documentation for arviz.from_dict
.
ArviZ.from_mcmcchains
— Functionfrom_mcmcchains(posterior::MCMCChains.Chains; kwargs...) -> InferenceData
from_mcmcchains(; kwargs...) -> InferenceData
from_mcmcchains(
posterior::MCMCChains.Chains,
posterior_predictive::Any,
predictions::Any,
log_likelihood::Any;
kwargs...
) -> InferenceData
Convert data in an MCMCChains.Chains
format into an InferenceData
.
Any keyword argument below without an an explicitly annotated type above is allowed, so long as it can be passed to convert_to_inference_data
.
Arguments
posterior::MCMCChains.Chains
: Draws from the posterior
Keywords
posterior_predictive::Any=nothing
: Draws from the posterior predictive distribution or name(s) of predictive variables inposterior
predictions::Any=nothing
: Out-of-sample predictions for the posterior.prior::Any=nothing
: Draws from the priorprior_predictive::Any=nothing
: Draws from the prior predictive distribution or name(s) of predictive variables inprior
observed_data::Dict{String,Array}=nothing
: Observed data on which theposterior
is conditional. It should only contain data which is modeled as a random variable. Keys are parameter names and values.constant_data::Dict{String,Array}=nothing
: Model constants, data included in the model which is not modeled as a random variable. Keys are parameter names and values.predictions_constant_data::Dict{String,Array}=nothing
: Constants relevant to the model predictions (i.e. newx
values in a linear regression).log_likelihood::Any=nothing
: Pointwise log-likelihood for the data. It is recommended to use this argument as a dictionary whose keys are observed variable names and whose values are log likelihood arrays.log_likelihood::String=nothing
: Name of variable inposterior
with log likelihoodslibrary=MCMCChains
: Name of library that generated the chainscoords::Dict{String,Vector}=Dict()
: Map from named dimension to named indicesdims::Dict{String,Vector{String}}=Dict()
: Map from variable name to names of its dimensionseltypes::Dict{String,DataType}=Dict()
: Apply eltypes to specific variables. This is used to assign discrete eltypes to discrete variables.
Returns
InferenceData
: The data with groups corresponding to the provided data
ArviZ.from_namedtuple
— Functionfrom_namedtuple(posterior::NamedTuple; kwargs...) -> InferenceData
from_namedtuple(posterior::Vector{<:NamedTuple}; kwargs...) -> InferenceData
from_namedtuple(posterior::Matrix{<:NamedTuple}; kwargs...) -> InferenceData
from_namedtuple(posterior::Vector{Vector{<:NamedTuple}}; kwargs...) -> InferenceData
from_namedtuple(
posterior::NamedTuple,
sample_stats::Any,
posterior_predictive::Any,
predictions::Any,
log_likelihood::Any;
kwargs...
) -> InferenceData
Convert a NamedTuple
or container of NamedTuple
s to an InferenceData
.
If containers are passed, they are flattened into a single NamedTuple
with array elements whose first dimensions correspond to the dimensions of the containers.
Arguments
posterior
: The data to be converted. It may be of the following types:::NamedTuple
: The keys are the variable names and the values are arrays with dimensions(nchains, ndraws, sizes...)
.::Vector{<:NamedTuple}
: Each element is aNamedTuple
from a chain withArray
/MonteCarloMeasurements.Particle
values with dimensions(ndraws, sizes...)
.::Matrix{<:NamedTuple}
: Each element is a single draw from a single chain, with array/scalar values with dimensionssizes
. The dimensions of the matrix container are(nchains, ndraws)
::Vector{Vector{<:NamedTuple}}
: The same as the above case.
Keywords
posterior_predictive::Any=nothing
: Draws from the posterior predictive distributionsample_stats::Any=nothing
: Statistics of the posterior sampling processpredictions::Any=nothing
: Out-of-sample predictions for the posterior.prior::Any=nothing
: Draws from the priorprior_predictive::Any=nothing
: Draws from the prior predictive distributionsample_stats_prior::Any=nothing
: Statistics of the prior sampling processobserved_data::Dict{String,Array}=nothing
: Observed data on which theposterior
is conditional. It should only contain data which is modeled as a random variable. Keys are parameter names and values.constant_data::Dict{String,Array}=nothing
: Model constants, data included in the model which is not modeled as a random variable. Keys are parameter names and values.predictions_constant_data::Dict{String,Array}=nothing
: Constants relevant to the model predictions (i.e. newx
values in a linear regression).log_likelihood::Any=nothing
: Pointwise log-likelihood for the data. It is recommended to use this argument as a dictionary whose keys are observed variable names and whose values are log likelihood arrays.library=nothing
: Name of library that generated the drawscoords::Dict{String,Vector}=Dict()
: Map from named dimension to named indicesdims::Dict{String,Vector{String}}=Dict()
: Map from variable name to names of its dimensions
Returns
InferenceData
: The data with groups corresponding to the provided data
Examples
using ArviZ
nchains, ndraws = 2, 10
data1 = (
x = rand(nchains, ndraws),
y = randn(nchains, ndraws, 2),
z = randn(nchains, ndraws, 3, 2),
)
idata1 = from_namedtuple(data1)
data2 = [(x = rand(ndraws), y = randn(ndraws, 2), z = randn(ndraws, 3, 2)) for _ = 1:nchains];
idata2 = from_namedtuple(data2)
data3 = [(x = rand(), y = randn(2), z = randn(3, 2)) for _ = 1:nchains, _ = 1:ndraws];
idata3 = from_namedtuple(data3)
data4 = [[(x = rand(), y = randn(2), z = randn(3, 2)) for _ = 1:ndraws] for _ = 1:nchains];
idata4 = from_namedtuple(data4)
ArviZ.from_netcdf
— FunctionSee documentation for arviz.from_netcdf
.
ArviZ.from_samplechains
— Functionfrom_samplechains(
posterior=nothing;
prior=nothing,
library=SampleChains,
kwargs...,
) -> InferenceData
Convert SampleChains samples to an InferenceData
.
Either posterior
or prior
may be a SampleChains.AbstractChain
or SampleChains.MultiChain
object.
For descriptions of remaining kwargs
, see from_namedtuple
.
ArviZ.hdi
— FunctionSee documentation for arviz.hdi
.
ArviZ.load_arviz_data
— FunctionSee documentation for arviz.load_arviz_data
.
ArviZ.loo
— FunctionSee documentation for arviz.loo
.
ArviZ.loo_pit
— FunctionSee documentation for arviz.loo_pit
.
ArviZ.make_ufunc
— FunctionSee documentation for arviz.make_ufunc
.
ArviZ.mcse
— FunctionSee documentation for arviz.mcse
.
ArviZ.plot_autocorr
— FunctionSee documentation for arviz.plot_autocorr
.
ArviZ.plot_bpv
— FunctionSee documentation for arviz.plot_bpv
.
ArviZ.plot_compare
— FunctionSee documentation for arviz.plot_compare
.
ArviZ.plot_density
— FunctionSee documentation for arviz.plot_density
.
ArviZ.plot_dist
— FunctionSee documentation for arviz.plot_dist
.
ArviZ.plot_dist_comparison
— FunctionSee documentation for arviz.plot_dist_comparison
.
ArviZ.plot_elpd
— FunctionSee documentation for arviz.plot_elpd
.
ArviZ.plot_energy
— FunctionSee documentation for arviz.plot_energy
.
ArviZ.plot_ess
— FunctionSee documentation for arviz.plot_ess
.
ArviZ.plot_forest
— FunctionSee documentation for arviz.plot_forest
.
ArviZ.plot_hdi
— FunctionSee documentation for arviz.plot_hdi
.
ArviZ.plot_kde
— FunctionSee documentation for arviz.plot_kde
.
ArviZ.plot_khat
— FunctionSee documentation for arviz.plot_khat
.
ArviZ.plot_loo_pit
— FunctionSee documentation for arviz.plot_loo_pit
.
ArviZ.plot_mcse
— FunctionSee documentation for arviz.plot_mcse
.
ArviZ.plot_pair
— FunctionSee documentation for arviz.plot_pair
.
ArviZ.plot_parallel
— FunctionSee documentation for arviz.plot_parallel
.
ArviZ.plot_posterior
— FunctionSee documentation for arviz.plot_posterior
.
ArviZ.plot_ppc
— FunctionSee documentation for arviz.plot_ppc
.
ArviZ.plot_rank
— FunctionSee documentation for arviz.plot_rank
.
ArviZ.plot_separation
— FunctionSee documentation for arviz.plot_separation
.
ArviZ.plot_trace
— FunctionSee documentation for arviz.plot_trace
.
ArviZ.plot_violin
— FunctionSee documentation for arviz.plot_violin
.
ArviZ.psislw
— Functionpsislw(log_weights, reff=1.0) -> (lw_out, kss)
Pareto smoothed importance sampling (PSIS).
This function is deprecated and is just a thin wrapper around psis
.
Arguments
log_weights
: Array of size(nobs, ndraws)
reff
: relative MCMC efficiency,ess / n
Returns
lw_out
: Smoothed log weightskss
: Pareto tail indices
ArviZ.r2_score
— FunctionSee documentation for arviz.r2_score
.
ArviZ.rhat
— FunctionSee documentation for arviz.rhat
.
ArviZ.to_netcdf
— FunctionSee documentation for arviz.to_netcdf
.
ArviZ.waic
— FunctionSee documentation for arviz.waic
.
ArviZ.with_interactive_backend
— Functionwith_interactive_backend(f; backend::Symbol = nothing)
Execute the thunk f
in a temporary interactive context with the chosen backend
, or provide no arguments to use a default.
Examples
idata = load_arviz_data("centered_eight")
plot_posterior(idata) # inline
with_interactive_backend() do
plot_density(idata) # interactive
end
plot_trace(idata) # inline
ArviZ.with_rc_context
— Functionwith_rc_context(f; rc = nothing, fname = nothing)
Execute the thunk f
within a context controlled by temporary rc params.
See rcParams
for supported params or to modify params long-term.
Examples
with_rc_context(fname = "pystan.rc") do
idata = load_arviz_data("radon")
plot_posterior(idata; var_names=["gamma"])
end
The plot would have settings from pystan.rc
.
A dictionary can also be passed to the context manager:
with_rc_context(rc = Dict("plot.max_subplots" => 1), fname = "pystan.rc") do
idata = load_arviz_data("radon")
plot_posterior(idata, var_names=["gamma"])
end
The rc
dictionary takes precedence over the settings loaded from fname
. Passing a dictionary only is also valid.
ArviZ.wrap_xarray_ufunc
— FunctionSee documentation for arviz.wrap_xarray_ufunc
.
StatsBase.summarystats
— Methodsummarystats(
data::InferenceData;
group = :posterior,
kwargs...,
) -> Union{Dataset,DataFrames.DataFrame}
summarystats(data::Dataset; kwargs...) -> Union{Dataset,DataFrames.DataFrame}
Compute summary statistics on data
.
Arguments
data::Union{Dataset,InferenceData}
: The data on which to compute summary statistics. Ifdata
is anInferenceData
, only the dataset corresponding togroup
is used.
Keywords
var_names::Vector{String}=nothing
: Names of variables to include in summaryinclude_circ::Bool=false
: Whether to include circular statisticsfmt::String="wide"
: Return format is eitherDataFrames.DataFrame
("wide", "long") orDataset
("xarray").round_to::Int=nothing
: Number of decimals used to round results. Usenothing
to return raw numbers.stat_funcs::Union{Dict{String,Function},Vector{Function}}=nothing
: A vector of functions or a dict of functions with function names as keys used to calculate statistics. By default, the mean, standard deviation, simulation standard error, and highest posterior density intervals are included. The functions will be given one argument, the samples for a variable as an array, The functions should operate on an array, returning a single number. For example,Statistics.mean
, orStatistics.var
would both work.extend::Bool=true
: Iftrue
, use the statistics returned bystat_funcs
in addition to, rather than in place of, the default statistics. This is only meaningful whenstat_funcs
is notnothing
.hdi_prob::Real=0.94
: HDI interval to compute. This is only meaningful whenstat_funcs
isnothing
.order::String="C"
: Iffmt
is "wide", use either "C" or "F" unpacking order.skipna::Bool=false
: Iftrue
, ignoresNaN
values when computing the summary statistics. It does not affect the behaviour of the functions passed tostat_funcs
.coords::Dict{String,Vector}=Dict()
: Coordinates specification to be used if thefmt
is"xarray"
.dims::Dict{String,Vector}=Dict()
: Dimensions specification for the variables to be used if thefmt
is"xarray"
.
Returns
Union{Dataset,DataFrames.DataFrame}
: Return type dicated byfmt
argument. Return value will contain summary statistics for each variable. Default statistics are:mean
sd
hdi_3%
hdi_97%
mcse_mean
mcse_sd
ess_bulk
ess_tail
r_hat
(only computed for traces with 2 or more chains)
Examples
using ArviZ
idata = load_arviz_data("centered_eight")
summarystats(idata; var_names=["mu", "tau"])
Other statistics can be calculated by passing a list of functions or a dictionary with key, function pairs:
using StatsBase, Statistics
function median_sd(x)
med = median(x)
sd = sqrt(mean((x .- med).^2))
return sd
end
func_dict = Dict(
"std" => x -> std(x; corrected = false),
"median_std" => median_sd,
"5%" => x -> percentile(x, 5),
"median" => median,
"95%" => x -> percentile(x, 95),
)
summarystats(idata; var_names = ["mu", "tau"], stat_funcs = func_dict, extend = false)
PSIS.PSISResult
— TypePSISResult
Result of Pareto-smoothed importance sampling (PSIS) using psis
.
Properties
log_weights
: un-normalized Pareto-smoothed log weightsweights
: normalized Pareto-smoothed weights (allocates a copy)pareto_shape
: Pareto $k=ξ$ shape parameternparams
: number of parameters inlog_weights
ndraws
: number of draws inlog_weights
nchains
: number of chains inlog_weights
reff
: the ratio of the effective sample size of the unsmoothed importance ratios and the actual sample size.ess
: estimated effective sample size of estimate of mean using smoothed importance samples (seeess_is
)log_weights_norm
: the logarithm of the normalization constant oflog_weights
tail_length
: length of the upper tail oflog_weights
that was smoothedtail_dist
: the generalized Pareto distribution that was fit to the tail oflog_weights
. Note that the tail weights are scaled to have a maximum of 1, sotail_dist * exp(maximum(log_ratios))
is the corresponding fit directly to the tail oflog_ratios
.
Diagnostic
The pareto_shape
parameter $k=ξ$ of the generalized Pareto distribution tail_dist
can be used to diagnose reliability and convergence of estimates using the importance weights [VehtariSimpson2021].
- if $k < \frac{1}{3}$, importance sampling is stable, and importance sampling (IS) and PSIS both are reliable.
- if $k ≤ \frac{1}{2}$, then the importance ratio distributon has finite variance, and the central limit theorem holds. As $k$ approaches the upper bound, IS becomes less reliable, while PSIS still works well but with a higher RMSE.
- if $\frac{1}{2} < k ≤ 0.7$, then the variance is infinite, and IS can behave quite poorly. However, PSIS works well in this regime.
- if $0.7 < k ≤ 1$, then it quickly becomes impractical to collect enough importance weights to reliably compute estimates, and importance sampling is not recommended.
- if $k > 1$, then neither the variance nor the mean of the raw importance ratios exists. The convergence rate is close to zero, and bias can be large with practical sample sizes.
See paretoshapeplot
for a diagnostic plot.
PSIS.ess_is
— Functioness_is(weights; reff=1)
Estimate effective sample size (ESS) for importance sampling over the sample dimensions.
Given normalized weights $w_{1:n}$, the ESS is estimated using the L2-norm of the weights:
\[\mathrm{ESS}(w_{1:n}) = \frac{r_{\mathrm{eff}}}{\sum_{i=1}^n w_i^2}\]
where $r_{\mathrm{eff}}$ is the relative efficiency of the log_weights
.
ess_is(result::PSISResult; bad_shape_missing=true)
Estimate ESS for Pareto-smoothed importance sampling.
ESS estimates for Pareto shape values $k > 0.7$, which are unreliable and misleadingly high, are set to missing
. To avoid this, set bad_shape_missing=false
.
PSIS.paretoshapeplot
— Functionparetoshapeplot(values; kwargs...)
paretoshapeplot!(values; kwargs...)
Plot shape parameters of fitted Pareto tail distributions for diagnosing convergence.
Arguments
values
: may be either a vector of Pareto shape parameters or aPSISResult
.
Keywords
showlines=false
: iftrue
, plot horizontal lines indicating relevant Pareto shape thresholds are drawn. SeePSISResult
for explanation of thresholds.backend::Symbol
: backend to use for plotting, defaulting to:Plots
, unless:Makie
is available.
All remaining keywords are passed to the plotting backend.
See psis
, PSISResult
.
Plots.jl or a Makie.jl backend must be loaded to use these functions.
Examples
using PSIS, Distributions, Plots
proposal = Normal()
target = TDist(7)
x = rand(proposal, 100, 1_000)
log_ratios = logpdf.(target, x) .- logpdf.(proposal, x)
result = psis(log_ratios)
Plot with Plots.jl.
using Plots
plot(result; showlines=true)
Plot with GLMakie.jl.
using GLMakie
plot(result; showlines=true)
PSIS.paretoshapeplot!
— Functionparetoshapeplot(values; kwargs...)
paretoshapeplot!(values; kwargs...)
Plot shape parameters of fitted Pareto tail distributions for diagnosing convergence.
Arguments
values
: may be either a vector of Pareto shape parameters or aPSISResult
.
Keywords
showlines=false
: iftrue
, plot horizontal lines indicating relevant Pareto shape thresholds are drawn. SeePSISResult
for explanation of thresholds.backend::Symbol
: backend to use for plotting, defaulting to:Plots
, unless:Makie
is available.
All remaining keywords are passed to the plotting backend.
See psis
, PSISResult
.
Plots.jl or a Makie.jl backend must be loaded to use these functions.
Examples
using PSIS, Distributions, Plots
proposal = Normal()
target = TDist(7)
x = rand(proposal, 100, 1_000)
log_ratios = logpdf.(target, x) .- logpdf.(proposal, x)
result = psis(log_ratios)
Plot with Plots.jl.
using Plots
plot(result; showlines=true)
Plot with GLMakie.jl.
using GLMakie
plot(result; showlines=true)
PSIS.psis
— Functionpsis(log_ratios, reff = 1.0; kwargs...) -> PSISResult
psis!(log_ratios, reff = 1.0; kwargs...) -> PSISResult
Compute Pareto smoothed importance sampling (PSIS) log weights [VehtariSimpson2021].
While psis
computes smoothed log weights out-of-place, psis!
smooths them in-place.
Arguments
log_ratios
: an array of logarithms of importance ratios, with one of the following sizes:(ndraws,)
: a vector of draws for a single parameter from a single chain(nparams, ndraws)
: a matrix of draws for a multiple parameter from a single chain(nparams, ndraws, nchains)
: an array of draws for multiple parameters from multiple chains, e.g. as might be generated with Markov chain Monte Carlo.
reff::Union{Real,AbstractVector}
: the ratio(s) of effective sample size oflog_ratios
and the actual sample sizereff = ess/(ndraws * nchains)
, used to account for autocorrelation, e.g. due to Markov chain Monte Carlo.
Keywords
improved=false
: Iftrue
, use the adaptive empirical prior of [Zhang2010]. Iffalse
, use the simpler prior of [ZhangStephens2009], which is also used in [VehtariSimpson2021].warn=true
: Iftrue
, warning messages are delivered
Returns
result
: aPSISResult
object containing the results of the Pareto-smoothing.
A warning is raised if the Pareto shape parameter $k ≥ 0.7$. See PSISResult
for details and paretoshapeplot
for a diagnostic plot.
PSIS.psis!
— Functionpsis(log_ratios, reff = 1.0; kwargs...) -> PSISResult
psis!(log_ratios, reff = 1.0; kwargs...) -> PSISResult
Compute Pareto smoothed importance sampling (PSIS) log weights [VehtariSimpson2021].
While psis
computes smoothed log weights out-of-place, psis!
smooths them in-place.
Arguments
log_ratios
: an array of logarithms of importance ratios, with one of the following sizes:(ndraws,)
: a vector of draws for a single parameter from a single chain(nparams, ndraws)
: a matrix of draws for a multiple parameter from a single chain(nparams, ndraws, nchains)
: an array of draws for multiple parameters from multiple chains, e.g. as might be generated with Markov chain Monte Carlo.
reff::Union{Real,AbstractVector}
: the ratio(s) of effective sample size oflog_ratios
and the actual sample sizereff = ess/(ndraws * nchains)
, used to account for autocorrelation, e.g. due to Markov chain Monte Carlo.
Keywords
improved=false
: Iftrue
, use the adaptive empirical prior of [Zhang2010]. Iffalse
, use the simpler prior of [ZhangStephens2009], which is also used in [VehtariSimpson2021].warn=true
: Iftrue
, warning messages are delivered
Returns
result
: aPSISResult
object containing the results of the Pareto-smoothing.
A warning is raised if the Pareto shape parameter $k ≥ 0.7$. See PSISResult
for details and paretoshapeplot
for a diagnostic plot.
Internal
ArviZ.BokehPlot
— TypeBokehPlot(::PyObject)
Loose wrapper around a Bokeh figure, mostly used for dispatch.
In most cases, use one of the plotting functions with backend=:bokeh
to create a BokehPlot
instead of using a constructor.
ArviZ.Dataset
— TypeDataset(::PyObject)
Dataset(; data_vars = nothing, coords = Dict(), attrs = Dict())
Loose wrapper around xarray.Dataset
, mostly used for dispatch.
Keywords
data_vars::Dict{String,Any}
: Dict mapping variable names toVector
: Data vector. Single dimension is named after variable.Tuple{String,Vector}
: Dimension name and data vector.Tuple{NTuple{N,String},Array{T,N}} where {N,T}
: Dimension names and data array.
coords::Dict{String,Any}
: Dict mapping dimension names to index names. Possible arguments has same form asdata_vars
.attrs::Dict{String,Any}
: Global attributes to save on this dataset.
In most cases, use convert_to_dataset
or convert_to_constant_dataset
or to create a Dataset
instead of directly using a constructor.
ArviZ.convert_arguments
— Methodconvert_arguments(f, args...; kwargs...) -> NTuple{2}
Convert arguments to the function f
before calling.
This function is used primarily for pre-processing arguments within macros before sending to arviz.
ArviZ.convert_result
— Methodconvert_result(f, result, args...)
Convert result of the function f
before returning.
This function is used primarily for post-processing outputs of arviz before returning. The args
are primarily used for dispatch.
ArviZ.convert_to_constant_dataset
— Functionconvert_to_constant_dataset(obj::Dict; kwargs...) -> Dataset
convert_to_constant_dataset(obj::NamedTuple; kwargs...) -> Dataset
Convert obj
into a Dataset
.
Unlike convert_to_dataset
, this is intended for containing constant parameters such as observed data and constant data, and the first two dimensions are not required to be the number of chains and draws.
Keywords
coords::Dict{String,Vector}
: Map from named dimension to index namesdims::Dict{String,Vector{String}}
: Map from variable name to names of its dimensionslibrary::Any
: A library associated with the data to add toattrs
.attrs::Dict{String,Any}
: Global attributes to save on this dataset.
ArviZ.convert_to_dataset
— Functionconvert_to_dataset(obj; group = :posterior, kwargs...) -> Dataset
Convert a supported object to a Dataset
.
In most cases, this function calls convert_to_inference_data
and returns the corresponding group
.
ArviZ.dataset_to_dict
— Functiondataset_to_dict(ds::Dataset) -> Tuple{Dict{String,Array},NamedTuple}
Convert a Dataset
to a dictionary of Array
s. The function also returns keyword arguments to dict_to_dataset
.
ArviZ.dict_to_dataset
— Functiondict_to_dataset(data::Dict{String,Array}; kwargs...) -> Dataset
Convert a dictionary with data and keys as variable names to a Dataset
.
Keywords
attrs::Dict{String,Any}
: Json serializable metadata to attach to the dataset, in addition to defaults.library::String
: Name of library used for performing inference. Will be attached to theattrs
metadata.coords::Dict{String,Array}
: Coordinates for the datasetdims::Dict{String,Vector{String}}
: Dimensions of each variable. The keys are variable names, values are vectors of coordinates.
Examples
using ArviZ
ArviZ.dict_to_dataset(Dict("x" => randn(4, 100), "y" => randn(4, 100)))
ArviZ.flatten
— Methodflatten(x)
If x
is an array of arrays, flatten into a single array whose dimensions are ordered with dimensions of the outermost container first and innermost container last.
ArviZ.from_cmdstanpy
— FunctionSee documentation for arviz.from_cmdstanpy
.
ArviZ.from_emcee
— FunctionSee documentation for arviz.from_emcee
.
ArviZ.from_numpyro
— FunctionSee documentation for arviz.from_numpyro
.
ArviZ.from_pymc3
— FunctionSee documentation for arviz.from_pymc3
.
ArviZ.from_pyro
— FunctionSee documentation for arviz.from_pyro
.
ArviZ.from_pystan
— FunctionSee documentation for arviz.from_pystan
.
ArviZ.from_tfp
— FunctionSee documentation for arviz.from_tfp
.
ArviZ.groupnames
— Methodgroupnames(data::InferenceData) -> Vector{Symbol}
Get the names of the groups (datasets) in data
.
ArviZ.groups
— Methodgroups(data::InferenceData) -> Dict{Symbol,Dataset}
Get the groups in data
as a dictionary mapping names to datasets.
ArviZ.namedtuple_of_arrays
— Methodnamedtuple_of_arrays(x::NamedTuple) -> NamedTuple
namedtuple_of_arrays(x::AbstractArray{NamedTuple}) -> NamedTuple
namedtuple_of_arrays(x::AbstractArray{AbstractArray{<:NamedTuple}}) -> NamedTuple
Given a container of NamedTuple
s, concatenate them, using the container dimensions as the dimensions of the resulting arrays.
Examples
using ArviZ
nchains, ndraws = 4, 100
data = [(x=rand(), y=randn(2), z=randn(2, 3)) for _ in 1:nchains, _ in 1:ndraws];
ntarray = ArviZ.namedtuple_of_arrays(data);
ArviZ.plot_ecdf
— FunctionSee documentation for arviz.plot_ecdf
.
ArviZ.styles
— Methodstyles() -> Vector{String}
Get all available matplotlib styles for use with use_style
ArviZ.summary
— Methodsummary(
data;
group = :posterior,
coords = nothing,
dims = nothing,
kwargs...,
) -> Union{Dataset,DataFrames.DataFrame}
Compute summary statistics on any object that can be passed to convert_to_dataset
.
Keywords
coords::Dict{String,Vector}=nothing
: Map from named dimension to named indices.dims::Dict{String,Vector{String}}=nothing
: Map from variable name to names of its dimensions.kwargs
: Keyword arguments passed tosummarystats
.
ArviZ.todataframes
— Methodtodataframes(df; index_name = nothing) -> DataFrames.DataFrame
Convert a Python pandas.DataFrame
or pandas.Series
into a DataFrames.DataFrame
.
If index_name
is not nothing
, the index is converted into a column with index_name
. Otherwise, it is discarded.
ArviZ.topandas
— Methodtopandas(::Type{:DataFrame}, df; index_name = nothing) -> PyObject
topandas(::Type{:Series}, df) -> PyObject
topandas(::Val{:ELPDData}, df) -> PyObject
Convert a DataFrames.DataFrame
to the specified pandas type.
If index_name
is not nothing
, the corresponding column is made the index of the returned dataframe.
ArviZ.use_style
— Methoduse_style(style::String)
use_style(style::Vector{String})
Use matplotlib style settings from a style specification style
.
The style name of "default" is reserved for reverting back to the default style settings.
ArviZ-specific styles are ["arviz-whitegrid", "arviz-darkgrid", "arviz-colors", "arviz-white"]
. To see all available style specifications, use styles()
.
If a Vector
of styles is provided, they are applied from first to last.
Base.write
— Methodwrite(io::IO, plot::BokehPlot)
write(filename::AbstractString, plot::BokehPlot)
Write the HTML representation of the Bokeh plot to the I/O stream or file.
ArviZ.@forwardfun
— Macro@forwardfun f
@forwardfun(f)
Wrap a function arviz.f
in f
, forwarding its docstrings.
Use convert_arguments
and convert_result
to customize what is passed to and returned from f
.
ArviZ.@forwardplotfun
— Macro@forwardplotfun f
@forwardplotfun(f)
Wrap a plotting function arviz.f
in f
, forwarding its docstrings.
This macro also ensures that outputs for the different backends are correctly handled. Use convert_arguments
and convert_result
to customize what is passed to and returned from f
.
- VehtariSimpson2021Vehtari A, Simpson D, Gelman A, Yao Y, Gabry J. (2021). Pareto smoothed importance sampling. arXiv:1507.02646v7 [stat.CO]
- VehtariSimpson2021Vehtari A, Simpson D, Gelman A, Yao Y, Gabry J. (2021). Pareto smoothed importance sampling. arXiv:1507.02646v7 [stat.CO]
- ZhangStephens2009Jin Zhang & Michael A. Stephens (2009) A New and Efficient Estimation Method for the Generalized Pareto Distribution, Technometrics, 51:3, 316-325, DOI: 10.1198/tech.2009.08017
- Zhang2010Jin Zhang (2010) Improving on Estimation for the Generalized Pareto Distribution, Technometrics, 52:3, 335-339, DOI: 10.1198/TECH.2010.09206
- VehtariSimpson2021Vehtari A, Simpson D, Gelman A, Yao Y, Gabry J. (2021). Pareto smoothed importance sampling. arXiv:1507.02646v7 [stat.CO]
- ZhangStephens2009Jin Zhang & Michael A. Stephens (2009) A New and Efficient Estimation Method for the Generalized Pareto Distribution, Technometrics, 51:3, 316-325, DOI: 10.1198/tech.2009.08017
- Zhang2010Jin Zhang (2010) Improving on Estimation for the Generalized Pareto Distribution, Technometrics, 52:3, 335-339, DOI: 10.1198/TECH.2010.09206