InformationGeometry.DataModel
— TypeIn addition to a DataSet
, a DataModel
contains the model as a function model(x,p)
and its derivative dmodel(x,p)
where x
denotes the x value of the data and p
is a vector of parameters on which the model depends. Crucially, dmodel
contains the derivatives of the model with respect to the parameters p
, not the x values.
InformationGeometry.DataSet
— TypeThe DataSet
type is a container for datapoints. It holds 3 vectors x, y, sigma
where the components of sigma
quantify the error bars associated with each measurement.
InformationGeometry.HyperCube
— TypeThe HyperCube
type has the fields vals::Vector{Vector}
, which stores the intervals which define the hypercube and dim::Int
, which gives the dimension. Overall it just offers a convenient and standardized way of passing domains for integration or plotting between functions without having to check that these domains are sensible every time. Examples for constructing HyperCube
s:
HyperCube([[1,3],[pi,2pi],[-500.0,100.0]])
HyperCube([[-1,1]])
HyperCube([-1,1])
HyperCube(LowerUpper([-1,-5],[0,-4]))
HyperCube(collect([-7,7.] for i in 1:3))
The HyperCube
type is closely related to the LowerUpper
type and they can be easily converted into each other. Examples for quantities that can be computed from and operations involving HyperCube
objects:
CubeVol(X)
TranslateCube(X,v::Vector)
CubeWidths(X)
InformationGeometry.LowerUpper
— TypeThe LowerUpper
type has a field L
and a field U
which respectively store the lower and upper boundaries of an N-dimensional Hypercube. It is very closely related to (and stores the same information as) the HyperCube
type. Examples for constructing LowerUpper
s:
LowerUpper([-1,-5,pi],[0,-4,2pi])
LowerUpper(HyperCube([[5,6],[-pi,0.5]]))
LowerUpper(collect(1:5),collect(15:20))
Examples for quantities that can be computed from and operations involving LowerUpper
objects:
CubeVol(X)
TranslateCube(X,v::Vector)
CubeWidths(X)
InformationGeometry.Plane
— TypeSpecifies a 2D plane in the so-called parameter form using 3 vectors.
InformationGeometry.AIC
— MethodAIC(DM::DataModel,p::Vector)
Calculates the Akaike Information Criterion given a parameter configuration p
.
InformationGeometry.ConfAlpha
— MethodConfAlpha(n::Real)
Probability volume outside of a confidence interval of level n⋅σ where σ is the standard deviation of a normal distribution.
InformationGeometry.ConfVol
— MethodConfVol(n::Real)
Probability volume contained in a confidence interval of level n⋅σ where σ is the standard deviation of a normal distribution.
InformationGeometry.CoverCubes
— MethodCoverCubes(A::HyperCube,B::HyperCube)
Return new HyperCube which covers two other given HyperCubes.
InformationGeometry.DataSpaceDist
— MethodDataSpaceDist(DM::DataModel,v::Vector) -> Real
Calculates the euclidean distance between a point v
in the data space and the data.
InformationGeometry.DistanceAlongGeodesic
— MethodDistanceAlongGeodesic(Metric::Function,sol::ODESolution,L::Real; tol=1e-14)
Calculates at which parameter value of the geodesic sol
the length L
is reached.
InformationGeometry.EvaluateEach
— MethodEvaluateEach(sols::Vector{Q}, Ts::Vector) where Q <: ODESolution
Evalues a family sols
of geodesics on a set of parameters Ts
. sols[1]
is evaluated at Ts[1]
, sols[2]
is evaluated at Ts[2]
and so on. The second half of the values respresenting the velocities is automatically truncated.
InformationGeometry.FindMLE
— FunctionFindMLE(DM::DataModel,start::Union{Bool,Vector}=false; Big::Bool=false)
Finds the maximum likelihood parameter configuration given a DataModel and optionally a starting configuration. Big=true
will return the value as a BigFloat
.
InformationGeometry.FisherMetric
— MethodFisherMetric(DM::DataModel, p::Vector{<:Real})
Computes the Fisher metric given a DataModel
and a parameter configuration p
under the assumption that the data is normally distributed around the true model. $F_{ab} = \int \mathrm{d}y \, p(y) \, (\partial_a \partial_b \ell)(x,\theta)$
InformationGeometry.GeodesicBetween
— MethodGeodesicBetween(Metric::Function,P::Vector{<:Real},Q::Vector{<:Real}; tol::Real=1e-10, meth=Tsit5())
Computes a geodesic between two given points on the parameter manifold and an expression for the metric.
InformationGeometry.GeodesicCrossing
— FunctionGeodesicCrossing(DM::DataModel,sol::ODESolution,Conf::Real=ConfVol(1); tol=1e-15)
Gives the parameter value of the geodesic sol
at which the confidence level Conf
is crossed.
InformationGeometry.GeodesicDistance
— MethodGeodesicDistance(Metric::Function,P::Vector{<:Real},Q::Vector{<:Real}; tol::Real=1e-10, meth=Tsit5())
Computes the length of a geodesic connecting the points P
and Q
.
InformationGeometry.GeodesicLength
— FunctionGeodesicLength(Metric::Function,sol::ODESolution, Endrange::Real=0.; fullSol::Bool=false, tol=1e-14)
Calculates the length of a geodesic sol
using the Metric
up to parameter value Endrange
.
InformationGeometry.KullbackLeibler
— FunctionKullbackLeibler(p::Function,q::Function,Domain::HyperCube=HyperCube([[-15,15]]); tol=1e-15, N::Int=Int(3e7), Carlo::Bool=(Domain.dim!=1))
Computes the Kullback-Leibler divergence between two probability distributions p
and q
over the Domain
. If Carlo=true
, this is done using a Monte Carlo Simulation with N
samples. If the Domain
is one-dimensional, the calculation is performed without Monte Carlo to a tolerance of ≈ tol
.
InformationGeometry.KullbackLeibler
— MethodKullbackLeibler(DM::DataModel,p::Vector,q::Vector)
Calculates KL divergence under assumption of normal error distribution around data.
InformationGeometry.NormalDist
— Methodfunction NormalDist(DM::DataModel,p::Vector) -> Distribution
Constructs either Normal
or MvNormal
type from Distributions.jl
using data and a parameter configuration. This makes the assumption, that the errors associated with the data are normal.
InformationGeometry.NumericChristoffel
— MethodNumericChristoffel(DM::DataModel, point::Vector; BigCalc::Bool=false)
NumericChristoffel(Metric::Function, point::Vector; BigCalc::Bool=false)
Calculates the Christoffel symbol at a point p
though finite differencing of the Metric
. Accurate to ≈ 3e-11. BigCalc=true
increases accuracy through BigFloat calculation.
InformationGeometry.OrthVF
— MethodOrthVF(DM::DataModel, p::Vector{<:Real}; Auto::Bool=false) -> Vector
Calculates a direction (in parameter space) in which the value of the log-likelihood does not change, given a parameter configuration p
. Auto=true
uses automatic differentiation to calculate the score.
InformationGeometry.OrthVF
— MethodOrthVF(DM::DataModel, PL::Plane, p::Vector{<:Real}; Auto::Bool=false) -> Vector
Calculates a direction (in parameter space) in which the value of the log-likelihood does not change, given a parameter configuration p
. If a Plane
is specified, the direction will be projected onto it. Auto=true
uses automatic differentiation to calculate the score.
InformationGeometry.PlotMatrix
— FunctionPlotMatrix(Mat::Matrix,MLE::Vector,N::Int=400)
Plots ellipse corresponding to a given covariance matrix which may additionally be offset by a vector MLE
.
InformationGeometry.PointwiseConfidenceBand
— MethodPointwiseConfidenceBand(DM::DataModel,sol::ODESolution,domain::HyperCube; N::Int=200)
Given a confidence interval sol
, the pointwise confidence band around the model prediction is computed for x values in domain
by evaluating the model on the boundary of the confidence interval.
InformationGeometry.Pullback
— MethodPullback(DM::DataModel, G::AbstractArray{<:Real,2}, point::Vector)
Pull-back of a (0,2)-tensor G
to the parameter manifold.
InformationGeometry.Pullback
— MethodPullback(DM::DataModel, omega::Vector{<:Real}, point::Vector) -> Vector
Pull-back of a covector to the parameter manifold.
InformationGeometry.Pushforward
— MethodPushforward(DM::DataModel, X::Vector, point::Vector)
Calculates the push-forward of a vector X
from the parameter manifold to the data space.
InformationGeometry.Ricci
— MethodRicci(Metric::Function, point::Vector; BigCalc::Bool=false)
Calculates the Ricci tensor by finite differencing of Metric
. BigCalc=true
increases accuracy through BigFloat calculation.
InformationGeometry.RicciScalar
— MethodRicciScalar(Metric::Function, point::Vector; BigCalc::Bool=false)
Calculates the Ricci Scalar by finite differencing of the Metric
. BigCalc=true
increases accuracy through BigFloat calculation.
InformationGeometry.Riemann
— MethodRiemann(Metric::Function, point::Vector; BigCalc::Bool=false)
Calculates the Riemann tensor by finite differencing of the Metric
. BigCalc=true
increases accuracy through BigFloat calculation.
InformationGeometry.Rsquared
— MethodRsquared(DM::DataModel,Fit::LsqFit.LsqFitResult) -> Real
Calculates the R² value of the fit result Fit
.
InformationGeometry.Score
— MethodScore(DM::DataModel,p::Vector{<:Real}; Auto::Bool=false)
Calculates the gradient of the log-likelihood with respect to a set of parameters p
. Auto=true
uses automatic differentiation.
InformationGeometry.ScoreDimN
— MethodCalculates the score of models with y vales of dim > 1.
InformationGeometry.Truncated
— MethodTruncated(sol::ODESolution) -> Function
Given a geodesic sol
, the second half of the components which represent the velocity are truncated off. The result gives the position of the geodesic as a function of the parameter. However, since it is no longer of type ODESolution
, one no longer has access to the fields sol.t
, sol.u
and so on.
InformationGeometry.VFRescale
— MethodVFRescale(ZeilenVecs::Array{<:Real,2},C::HyperCube;scaling=0.85)
Rescale vector to look good in 2D plot.
InformationGeometry.WilksTest
— FunctionWilksTest(DM::DataModel, p::Vector{<:Real},MLE::Vector{<:Real},ConfVol=ConfVol(1)) -> Bool
Checks whether a given parameter configuration p
is within a confidence interval of level ConfVol
using Wilks' theorem. This makes the assumption, that the likelihood has the form of a normal distribution, which is asymptotically correct in the limit that the number of datapoints is infinite.
InformationGeometry.WilksTestPrepared
— FunctionWilksTestPrepared(DM::DataModel, p::Vector{<:Real}, LogLikeMLE::Real, ConfVol=ConfVol(1)) -> Bool
Checks whether a given parameter configuration p
is within a confidence interval of level ConfVol
using Wilks' theorem. This makes the assumption, that the likelihood has the form of a normal distribution, which is asymptotically correct in the limit that the number of datapoints is infinite. To simplify the calculation, LogLikeMLE
accepts the value of the log-likelihood evaluated at the MLE.
InformationGeometry.likelihood
— Methodlikelihood(DM::DataModel,p::Vector)
Calculates the Likelihood given a DataModel and a parameter configuration p
.
StatsBase.loglikelihood
— Methodloglikelihood(DM::DataModel, p::Vector)
Calculates the logarithm of the Likelihood given a DataModel and a parameter configuration p
.