Probability & Statistics

High Level Functions

Probability, NProbability

Expectation, NExpectation

RandomVariate

EstimatedDistribution, DistributionFitTest

MomentConvert, MomentEvaluate

Distribution Properties

PDF, CDF, SurvivalFunction (SF), HazardFunction (HF), Quantile (QF), InverseCDF, InverseSurvivalFunction

Moment / MGF, CentralMoment / CMGF, FactorialMoment / FMGF, Cumulant / CGF, CF

Mean, Variance, StandardDeviation, Skewness, Kurtosis, Covariance, Correlation, ...

Median, Quartiles, QuartileDeviation, QuartileSkewness, ...

Likelihood, LogLikelihood, ...

Parametric Distributions

Distributions constructed from parameters.

Univariate Discrete

Univariate Continuous

Multivariate Discrete

Multivariate Continuous

Nonparametric Distributions

Distributions constructed “directly” from data.

EmpiricalDistribution

HistogramDistribution

SmoothKernelDistribution

SurvivalDistribution

Formula Distributions

Distributions constructed from formulas.

ProbabilityDistribution (from PDF, CDF, SF)

Derived Distributions

Distributions constructed from other distributions.

TransformedDistribution

OrderDistribution

MixtureDistribution

ParameterMixtureDistribution

CensoredDistribution

TruncatedDistribution

MarginalDistribution

ProductDistribution

CopulaDistribution

Statistical Visualization

Understand and compare distributions.

Plot, Plot3D, DiscretePlot, DiscretePlot3D, ...

Histogram, Histogram3D, DensityHistogram, PairedHistogram, ...

SmoothHistogram, SmoothHistogram3D, SmoothDensityHistogram, ...

QuantilePlot, ProbabilityPlot, ProbabilityScalePlot,  ...

BoxWhiskerChart, DistributionChart, ...

Statistical Data

Access to data sources

Curated data (FinancialData, WeatherData, CountryData, ...)

Statistical data (ExampleData)

Import formats (about 160 formats)

Database access (MySQL, Oracle, ...)

Sommario

Probability_

FramedStatement[Grid[{{HoldForm[Probability[a,x\[Distributed]𝒟]],HoldForm[NProbability[a,x\[Distributed]𝒟]]},{HoldForm[Probability[a\[Conditioned]b,x\[Distributed]𝒟]],HoldForm[NProbability[a\[Conditioned]b,x\[Distributed]𝒟]]}},Alignment->{Left}]]

Definitions

FormattedDefinitionGrid[{{HoldForm@((∑_(ξ=-∞)^∞)Boole[a]f(ξ) ),HoldForm@((∑_((Subscript[ξ, 1])=-∞)^∞)(∑_((Subscript[ξ, n])=-∞)^∞)Boole[a] f((Subscript[ξ, 1]),…,(Subscript[ξ, n])))},{HoldForm@((∫_(-∞)^∞)Boole[a]f(ξ) \[DifferentialD]ξ),HoldForm@((∫_(-∞)^∞)(∫_(-∞)^∞)Boole[a] f((Subscript[ξ, 1]),…,(Subscript[ξ, n]))\[DifferentialD](Subscript[ξ, n])\[DifferentialD](Subscript[ξ, 1]))}}]

Basics

Discrete univariate probabilities:

Probability[x≥1,x\[Distributed]BinomialDistribution[n,p]]

Nonlinear conditions:

Probability[(((x-3))^2)≥x,x\[Distributed]ZipfDistribution[ρ]]

Logical combinations:

Probability[k^2>10||k<2,k\[Distributed]PoissonDistribution[λ]]

Multivariate probabilities:

Probability[X≤1/2&&Y≤1,{X,Y}\[Distributed]DirichletDistribution[{1,2,3}]]

Conditional probabilities:

Probability[(x^2)>12\[Conditioned]x>2,x\[Distributed]PoissonDistribution[λ]]

Applications

An airplane engine fails with probability p. Compute the probability that no more then 2 out of 4 fail:

prob1=Probability[x≤2,x\[Distributed]BinomialDistribution[4,p]]

Compute the probability that no more then 1 out of 2 fail:

prob2=Probability[x≤1,x\[Distributed]BinomialDistribution[2,p]]

When are four engines better then two engines?

Reduce[prob1>prob2&&0<p<1,p]

Plot[{Tooltip[prob1,"4 engines"],Tooltip[prob2,"2 engines"]},{p,0,1}]

A company manufactures nails with length normally distributed and a mean of 0.5 inches. If  50% of the produced nails have lengths between 0.495 and 0.505, find the standard deviation:

prob=Probability[0.495<x<0.505,x\[Distributed]NormalDistribution[0.5,σ]]

LogLinearPlot[prob,{σ,10^-3,1},GridLines->Automatic,Frame->True]

NSolve[prob==0.5&&10^-3<σ<10^-1,σ]

Expectation

FramedStatement[Grid[{{HoldForm[Expectation[g[x],x\[Distributed]𝒟]],HoldForm[NExpectation[g[x],x\[Distributed]𝒟]]},{HoldForm[Expectation[g[x]\[Conditioned]b,x\[Distributed]𝒟]],HoldForm[NExpectation[g[x]\[Conditioned]b,x\[Distributed]𝒟]]}},Alignment->{Left}]]

Definitions

FormattedDefinitionGrid[{{HoldForm@((∑_(ξ=-∞)^∞)g[ξ]f(ξ) ),HoldForm@((∑_((Subscript[ξ, 1])=-∞)^∞)(∑_((Subscript[ξ, n])=-∞)^∞)g[(Subscript[ξ, 1]),…,(Subscript[ξ, n])] f((Subscript[ξ, 1]),…,(Subscript[ξ, n])))},{HoldForm@((∫_(-∞)^∞)g[ξ]f(ξ) \[DifferentialD]ξ),HoldForm@((∫_(-∞)^∞)(∫_(-∞)^∞)g[(Subscript[ξ, 1]),…,(Subscript[ξ, n])] f((Subscript[ξ, 1]),…,(Subscript[ξ, n]))\[DifferentialD](Subscript[ξ, n])\[DifferentialD](Subscript[ξ, 1]))}}]

Basics

Standard polynomial moments:

Expectation[(k^2),k\[Distributed]PoissonDistribution[μ]]

General expectations:

Expectation[(2^(-k)),k\[Distributed]PoissonDistribution[μ]]

Sometimes only special moments exists:

Expectation[Exp[-Abs[x]],x\[Distributed]CauchyDistribution[1,2]]

N[%]

NExpectation[Exp[-Abs[x]],x\[Distributed]CauchyDistribution[1,2]]

Application

A basketball player shoots free throws until he hits 4 of them. His probability of scoring in any one of them is 0.7. Find the number of shots the player is expected to shoot:

Expectation[x+4,x\[Distributed]NegativeBinomialDistribution[4,0.7]]

Expectation[x,x\[Distributed]PascalDistribution[4,0.7]]

Assume the probability of fouling for each minute interval is 0.1 independently. A basketball player fouls out after 6 fouls. Find expected playing time until foul out:

Expectation[x+6,x\[Distributed]NegativeBinomialDistribution[6,0.1]]

Distribution Properties_

Distribution Functions

FormulaGallery[{{HoldForm[PDF[𝒟,x]]},{HoldForm[CDF[𝒟,x]]},{HoldForm[SurvivalFunction[𝒟,x]]},{HoldForm[HazardFunction[𝒟,x]]},{HoldForm[Quantile[𝒟,q]]},{HoldForm[InverseCDF[𝒟,q]]},{HoldForm[InverseSurvivalFunction[𝒟,q]]}}]

Moments

FormulaGallery[Transpose@{{HoldForm[Mean[𝒟]],HoldForm[Variance[𝒟]],HoldForm[StandardDeviation[𝒟]],HoldForm[Skewness[𝒟]],HoldForm[Kurtosis[𝒟]],HoldForm[Covariance[𝒟]],HoldForm[Correlation[𝒟]],HoldForm[Moment[𝒟,r]],HoldForm[CentralMoment[𝒟,r]],HoldForm[Cumulant[𝒟,r]],HoldForm[FactorialMoment[𝒟,r]]}}]

Quantiles

FormulaGallery[Transpose@{{HoldForm[Median[𝒟]],HoldForm[InterquartileRange[𝒟]],HoldForm[QuartileDeviation[𝒟]],HoldForm[Quartiles[𝒟]]}}]

Generating Functions

FormulaGallery[Transpose@{{HoldForm[MomentGeneratingFunction[𝒟,t]],HoldForm[CentralMomentGeneratingFunction[𝒟,t]],
HoldForm[CumulantGeneratingFunction[𝒟,t]],HoldForm[FactorialMomentGeneratingFunction[𝒟,t]],
HoldForm[CharacteristicFunction[𝒟,t]]}}]

Simulation, Estimation & Testing

FormulaGallery[Transpose@{{HoldForm[Likelihood[𝒟,data]],HoldForm[LogLikelihood[𝒟,data]],HoldForm[EstimatedDistribution[data,𝒟]],HoldForm[FindDistributionParameters[data,𝒟]],HoldForm[DistributionFitTest[data,𝒟]],HoldForm[QuantilePlot[data,𝒟]],HoldForm[ProbabilityPlot[data,𝒟]],HoldForm[ProbabilityScalePlot[data,𝒟]],HoldForm[RandomVariate[𝒟,n]]}}]

Distribution Functions

DistributionPlot[NormalDistribution[],{-3,3}]//CellDisplay

Definitions

CellDisplay@Style[TabView[{"PDF"->DefinitionGrid[{{HoldForm@Probability[X==k,…],HoldForm@Probability[(Subscript[X, 1])==(Subscript[k, 1])∧⋯∧(Subscript[X, n])==(Subscript[k, n]),…]},{HoldForm@Limit[((Probability[x-ϵ<X≤x,…])/ϵ),ϵ->0],HoldForm@Limit[((Probability[(Subscript[x, 1])-(Subscript[ϵ, 1])<(Subscript[X, 1])≤(Subscript[x, 1])∧⋯∧(Subscript[x, n])-(Subscript[ϵ, n])<(Subscript[X, n])≤(Subscript[x, n]),…])/((Subscript[ϵ, 1])⋯ (Subscript[ϵ, n]))),(Subscript[ϵ, 1])->0,…,(Subscript[ϵ, n])->0]}}],
"CDF"->DefinitionGrid[{{HoldForm@Probability[X≤k,…],HoldForm@Probability[(Subscript[X, 1])≤(Subscript[k, 1])∧⋯∧(Subscript[X, n])≤(Subscript[k, n]),…]},{HoldForm@Probability[X≤x,…],HoldForm@Probability[(Subscript[X, 1])≤(Subscript[x, 1])∧⋯∧(Subscript[X, n])≤(Subscript[x, n]),…]}}],
"SF"->DefinitionGrid[{{HoldForm@Probability[X>k,…],HoldForm@Probability[(Subscript[X, 1])>(Subscript[k, 1])∧⋯∧(Subscript[X, n])>(Subscript[k, n]),…]},{HoldForm@Probability[X>x,…],HoldForm@Probability[(Subscript[X, 1])>(Subscript[x, 1])∧⋯∧(Subscript[X, n])>(Subscript[x, n]),…]}}],
"HF"->DefinitionGrid[{{HoldForm@Probability[X==k\[Conditioned]X≥k,…],HoldForm@Probability[((Subscript[X, 1])==(Subscript[k, 1])∧⋯∧(Subscript[X, n])==(Subscript[k, n]))\[Conditioned]((Subscript[X, 1])≥(Subscript[k, 1])∧⋯∧(Subscript[X, n])≥(Subscript[k, n])),…]},{HoldForm@Limit[((Probability[x<X≤x+ϵ\[Conditioned]X≥x,…])/ϵ),ϵ->0],HoldForm@Limit[((Probability[((Subscript[x, 1])<(Subscript[X, 1])≤(Subscript[x, 1])+(Subscript[ϵ, 1])∧⋯∧(Subscript[x, n])<(Subscript[X, n])≤(Subscript[x, n])+(Subscript[ϵ, n]))\[Conditioned]((Subscript[X, 1])≥(Subscript[x, 1])∧⋯∧(Subscript[X, n])≥(Subscript[x, n])),…])/((Subscript[ϵ, 1])⋯ (Subscript[ϵ, n]))),(Subscript[ϵ, 1])->0,…,(Subscript[ϵ, n])->0]}}]}],"Text"]

Basics

Continuous univariate:

DistributionPlot[GompertzMakehamDistribution[1,2],{0,2}]

Discrete univariate:

DistributionPlot[BinomialDistribution[10,2/3],{0,11}]

Continuous bivariate:

DistributionPlot3D[BinormalDistribution[1/2],{-3,3},{-3,3}]

Discrete bivariate:

DistributionPlot3D[MultivariatePoissonDistribution[2,{1,2}],{0,10},{0,10}]

These can be accessed individually for each distribution:

PDF[ExponentialDistribution[λ],x]

CDF[ExponentialDistribution[λ],x]

SurvivalFunction[ExponentialDistribution[λ],x]

HazardFunction[ExponentialDistribution[λ],x]

Moments

FramedStatement[Column[{HoldForm[Moment[𝒟,k]],HoldForm[Cumulant[𝒟,k]],HoldForm[CentralMoment[𝒟,k]],HoldForm[FactorialMoment[𝒟,k]]}]]

Definitions

Moments and generating functions ... starting from moments:

CellDisplay@TabView[{"Moment"->MomentDefinitionGrid[{{HoldForm[(Subscript[μ, n])==Expectation[(X^n),…]],HoldForm[MGF[t]==GeneratingFunction[((Subscript[μ, n])/(n!)),n,t]]},{HoldForm[(Subscript[μ, (Subscript[n, 1]),…,(Subscript[n, k])])==Expectation[((Subscript[X, 1])^(Subscript[n, 1]))⋯ ((Subscript[X, k])^(Subscript[n, k])),…]],HoldForm[MGF[(Subscript[t, 1]),…,(Subscript[t, k])]==GeneratingFunction[((Subscript[μ, (Subscript[n, 1]),…,(Subscript[n, k])])/((Subscript[n, 1])!⋯ (Subscript[n, k])!)),{(Subscript[n, 1]),…,(Subscript[n, k])},{(Subscript[t, 1]),…,(Subscript[t, k])}]]}}],"CentralMoment"->MomentDefinitionGrid[{{HoldForm[(Subscript[(OverTilde[μ]), n])==Expectation[(((X-μ))^n),…]],HoldForm[CMGF[t]==GeneratingFunction[((Subscript[(OverTilde[μ]), n])/(n!)),n,t]]},{HoldForm[(Subscript[(OverTilde[μ]), (Subscript[n, 1]),…,(Subscript[n, k])])==Expectation[((((Subscript[X, 1])-(Subscript[μ, (Subscript[X, 1])])))^(Subscript[n, 1]))⋯ ((((Subscript[X, k])-(Subscript[μ, (Subscript[X, k])])))^(Subscript[n, k])),…]],HoldForm[CMGF[(Subscript[t, 1]),…,(Subscript[t, k])]==GeneratingFunction[((Subscript[(OverTilde[μ]), (Subscript[n, 1]),…,(Subscript[n, k])])/((Subscript[n, 1])!⋯ (Subscript[n, k])!)),{(Subscript[n, 1]),…,(Subscript[n, k])},{(Subscript[t, 1]),…,(Subscript[t, k])}]]}}],"FactorialMoment"->MomentDefinitionGrid[{{HoldForm[(Subscript[(OverBar[μ]), n])==Expectation[FactorialPower[X,n],…]],HoldForm[FMGF[t]==GeneratingFunction[((Subscript[(OverBar[μ]), n])/(n!)),n,t-1]]},{HoldForm[(Subscript[(OverBar[μ]), (Subscript[n, 1]),…,(Subscript[n, k])])==Expectation[FactorialPower[(Subscript[X, 1]),(Subscript[n, 1])]⋯ FactorialPower[(Subscript[X, k]),(Subscript[n, k])],…]],HoldForm[FMGF[(Subscript[t, 1]),…,(Subscript[t, k])]==GeneratingFunction[((Subscript[(OverBar[μ]), (Subscript[n, 1]),…,(Subscript[n, k])])/((Subscript[n, 1])!⋯ (Subscript[n, k])!)),{(Subscript[n, 1]),…,(Subscript[n, k])},{(Subscript[t, 1])-1,…,(Subscript[t, k])-1}]]}}]}]

Generating functions and moments ... starting from generating functions:

CellDisplay@Style[TabView[{"MGF"->MomentDefinitionGrid[{{HoldForm[(Subscript[μ, n])==Derivative[n][MGF][0]],HoldForm[MGF[t]==Expectation[Exp[t x],…]]},{HoldForm[(Subscript[μ, (Subscript[n, 1]),…,(Subscript[n, k])])==Derivative[(Subscript[n, 1]),…,(Subscript[n, k])][MGF][0,…,0]],HoldForm[MGF[(Subscript[t, 1]),…,(Subscript[t, k])]==Expectation[Exp[(Subscript[t, 1])(Subscript[x, 1])+⋯+(Subscript[t, k])(Subscript[x, k])],…]]}}],"CMGF"->MomentDefinitionGrid[{{HoldForm[(Subscript[(OverTilde[μ]), n])==Derivative[n][CMGF][0]],HoldForm[CMGF[t]==Expectation[Exp[t (x-μ)],…]]},{HoldForm[(Subscript[(OverTilde[μ]), (Subscript[n, 1]),…,(Subscript[n, k])])==Derivative[(Subscript[n, 1]),…,(Subscript[n, k])][CMGF][0,…,0]],HoldForm[CMGF[(Subscript[t, 1]),…,(Subscript[t, k])]==Expectation[Exp[(Subscript[t, 1])((Subscript[x, 1])-(Subscript[μ, (Subscript[x, 1])]))+⋯+(Subscript[t, k])((Subscript[x, k])-(Subscript[μ, (Subscript[x, k])]))],…]]}}],"FMGF"->MomentDefinitionGrid[{{HoldForm[(Subscript[(OverBar[μ]), n])==Derivative[n][FMGF][1]],HoldForm[FMGF[t]==Expectation[(t^x),…]]},{HoldForm[(Subscript[(OverBar[μ]), (Subscript[n, 1]),…,(Subscript[n, k])])==Derivative[(Subscript[n, 1]),…,(Subscript[n, k])][FMGF][1,…,1]],HoldForm[FMGF[(Subscript[t, 1]),…,(Subscript[t, k])]==Expectation[((Subscript[t, 1])^(Subscript[x, 1]))⋯ ((Subscript[t, k])^(Subscript[x, k])),…]]}}],
"CGF"->MomentDefinitionGrid[{{HoldForm[(Subscript[κ, n])==Derivative[n][CGF][0]],HoldForm[CGF[t]==Log[MGF[t]]]},{HoldForm[(Subscript[κ, (Subscript[n, 1]),…,(Subscript[n, k])])==Derivative[(Subscript[n, 1]),…,(Subscript[n, k])][CGF][0,…,0]],HoldForm[CGF[(Subscript[t, 1]),…,(Subscript[t, k])]==Log[MGF[(Subscript[t, 1]),…,(Subscript[t, k])]]]}}],
"CF"->MomentDefinitionGrid[{{…,HoldForm[CF[t]==Expectation[Exp[I t x],…]]},{…,HoldForm[CF[(Subscript[t, 1]),…,(Subscript[t, k])]==Expectation[Exp[I ((Subscript[t, 1])(Subscript[x, 1])+⋯+(Subscript[t, k])(Subscript[x, k]))],…]]}}]}],"Text"]

Basics

Moments for distributions:

Table[Moment[ExponentialDistribution[λ],k],{k,5}]

MomentTabView[NormalDistribution[μ,σ],5]

Moment conversions:

MomentConvert[CentralMoment[2],Moment]//TraditionalForm

CellDisplay@Manipulate[TraditionalForm[f[n]]->TraditionalForm@MomentConvert[f[n],g],{{f,Moment,"source"},{Moment,CentralMoment,FactorialMoment,Cumulant}},{{g,CentralMoment,"target"},{Moment,CentralMoment,FactorialMoment,Cumulant}},{{n,2},1,20,1}]

Applications

Build a method of moments estimator:

{(Subscript[m, 1])==Moment[NormalDistribution[μ,σ],1],(Subscript[m, 2])==Simplify@Moment[NormalDistribution[μ,σ],2]}

Solve[%,{μ,σ}]

data=RandomVariate[NormalDistribution[1,2],(10^4)];

NSolve[{Moment[data,1]==Moment[NormalDistribution[μ,σ],1],Moment[data,2]==Moment[NormalDistribution[μ,σ],2],σ>0},{μ,σ}]

Guess the closed form for moments:

Table[Moment[ExponentialDistribution[λ],k],{k,5}]

FindSequenceFunction[%,k]

Moment[ExponentialDistribution[λ],k]

FullSimplify[%%-%]

Find the moments from an MGF:

MomentGeneratingFunction[ExponentialDistribution[λ],t]

Table[Derivative[n][Function[t,(λ/(-t+λ))]][0],{n,5}]

SeriesCoefficient[n!(λ/(-t+λ)),{t,0,n}]

Simulation_

(*
data=RandomNumber[NormalDistribution[0,1],10^3];
datap=ListLinePlot[data,Ticks->None,ImageSize->Tiny];
distp=Plot[{0,PDF[NormalDistribution[0,1],x]},{x,-3,3},Filling->Axis,Ticks->None,ImageSize->Tiny];
*)
FramedStatement[HoldForm[RandomVariate[,…]⟹ ]]

Basics

Given a distribution we can simulate it:

data=RandomVariate[GumbelDistribution[1,3],10^4];

Show[
Histogram[data,25,"PDF"],
Plot[PDF[GumbelDistribution[1,3],x],{x,-15,8},PlotStyle->Thick]]

Similarly for bivariate random variates:

data=RandomVariate[BinormalDistribution[{1,2},{1.5,2},0.6],10^4];

Show[Histogram3D[data,{{-5,6,1},{-3,7,1}},"PDF",ChartBaseStyle->Opacity[0.5]],
Plot3D[PDF[BinormalDistribution[{1,2},{1.5,2},0.6],{x,y}],{x,-5,6},{y,-3,7},PlotStyle->Hue[.15,.7,.8],Mesh->None]]

Another example

dist=TruncatedDistribution[{1,4},TriangularDistribution[{0,7}]];

data=RandomVariate[dist,10^5];

Show[
Histogram[data,20,"ProbabilityDensity",PlotRange-> All],
Plot[PDF[dist,x],{x,0,5},PlotStyle->{Thickness[0.01],Red},Exclusions->None],AxesOrigin->{0,0}]

Mean, variance, and kurtosis:

{N@Mean[dist],Mean[data]}

{N@Variance[dist],Variance[data]}

{N@Kurtosis[dist],Kurtosis[data]}

Applications

The average number of accidents in a city is 100 per day. Simulate the daily accidents:

RandomVariate[PoissonDistribution[100],30]

ListPlot[%,Filling->Axis]

Estimation_

(*
data=RandomNumber[NormalDistribution[0,1],10^3];
datap=ListLinePlot[data,Ticks->None,ImageSize->Tiny];
distp=Plot[{0,PDF[NormalDistribution[0,1],x]},{x,-3,3},Filling->Axis,Ticks->None,ImageSize->Tiny];
distfamp=Plot[Evaluate[Table[PDF[NormalDistribution@@θ,x],{θ,{{-1,1/2},{0,1},{1,2/3},{1/2,2},{0,1/3},{2/3,1}}}]],{x,-3,3},Filling->Axis,Ticks->None,PlotRange->All,ImageSize->Tiny];
*)


FramedStatement[HoldForm[EstimatedDistribution[,]⟹]]

Basics

And given outcomes we can estimate the distribution:

data=RandomVariate[PoissonDistribution[100],10^3];

dist=EstimatedDistribution[data,PoissonDistribution[μ]]

Compare data with estimated distribution:

QuantilePlot[data,dist]

Applications

Get the distributions word length for some languages:

languages={"Dutch","English","French","German","Hebrew","Polish","Russian","Spanish","Swedish"};

worddata=Table[StringLength/@DictionaryLookup[{l,All}],{l,languages}];

Estimate binomial models:

binom=ParallelTable[EstimatedDistribution[i,BinomialDistribution[n,p]],{i,worddata}]//Quiet

Compare results:

comp=ParallelTable[Show[Histogram[worddata[[i]],{Range[25]-1/2},"PDF",PlotLabel->languages[[i]]],DiscretePlot[PDF[binom[[i]],x],{x,0,25},PlotRange->All,PlotStyle->PointSize[.02]]],{i,Length[languages]}];

GraphicsGrid[Partition[comp,3],ImageSize->Large]

Hypothesis Testing

(*
data=RandomNumber[NormalDistribution[0,1],10^3];
datap=ListLinePlot[data,Ticks->None,ImageSize->Tiny];
distp=Plot[{0,PDF[NormalDistribution[0,1],x]},{x,-3,3},Filling->Axis,Ticks->None,ImageSize->Tiny];
distfamp=Plot[Evaluate[Table[PDF[NormalDistribution@@θ,x],{θ,{{-1,1/2},{0,1},{1,2/3},{1/2,2},{0,1/3},{2/3,1}}}]],{x,-3,3},Filling->Axis,Ticks->None,PlotRange->All,ImageSize->Tiny];
*)


FramedStatement[HoldForm[DistributionFitTest[,]⟹]]

Basics

Get the probability that the data came from a normal distribution:

data=RandomVariate[NormalDistribution[],10^4];

DistributionFitTest[data]

ProgressIndicator[%]

Get a report from different goodness-of-fit tests:

ℋ=DistributionFitTest[data,Automatic,"HypothesisTestData"];

ℋ["TestDataTable",All]

BarChart[ℋ["TestData",All][[All,2]],ChartLabels->Placed[ℋ["AllTests"],Center],BarOrigin->Left]

The probability is itself random and UniformDistribution[{0,1}]:

data=RandomVariate[𝒹=NormalDistribution[],{1000,50}];

pdata=DistributionFitTest/@data;

Histogram[pdata,Automatic,"PDF"]

Parametric Distributions

Distributions constructed from parameters: "Parametric Distributions Palette"

Coverage

Compared to some close competitors:

Comparing both properties and distributions:

Nonparametric Distributions

FramedLargeStatement@Grid[{
{EmpiricalDistribution,"distribution from all data"},
{HistogramDistribution,"distribution from binning data"},
{SmoothKernelDistribution,"distribution from kernel smoothed data"},
{SurvivalDistribution,"empirical distribution for censored data"}},
Alignment->{Left,Baseline}]

Nonparametric Distributions – HistogramDistribution

(* 
𝒟ref=ComponentMixtureDistribution[{3/4,1/4},{NormalDistribution[],NormalDistribution[3/2,1/3]}];
datap=ListLinePlot[RandomNumber[𝒟ref,10^3],Ticks->None,ImageSize->Tiny];
𝒟dist=Quiet@HistogramDistribution[RandomNumber[𝒟ref,10^4],{-3,3,0.1}];
distp=Plot[{0,PDF[𝒟dist,x]},{x,-3,3},Filling->Axis,PlotPoints->200,Ticks->None,ImageSize->Tiny]; 
*)



FramedStatement[HoldForm[HistogramDistribution[,…]⟹ ]]

(*
𝒟ref=ComponentMixtureDistribution[{1/3,2/3},{BinormalDistribution[{-1,-1},{1,1},1/2],BinormalDistribution[{1,1},{1,1},-1/2]}];
sdata=Graphics[{ColorData[1,1],Point@RandomNumber[𝒟ref,10^3]},Ticks->None,ImageSize->Tiny]
data=RandomNumber[𝒟ref,10^4];
𝒟dist=Quiet@HistogramDistribution[RandomNumber[𝒟ref,10^4],{{-3,3,0.5},{-3,3,0.5}}];
distp=Histogram3D[{{{0,0},{0,1}},data},{{-3,3,0.5},{-3,3,0.5}},Ticks->None]
*)


FramedStatement[HoldForm[HistogramDistribution[,…]⟹ ]]

Basics 1D

Normally we think of histogram as a graphical technique for approximately visualizing distributions:

Histogram[RandomVariate[NormalDistribution[],10^4]]

But we can just as well think of it as an actual distribution:

𝒟=HistogramDistribution[RandomVariate[NormalDistribution[],10^4],6]

Behaves just like any other distribution:

PDF[𝒟,x]

CDF[𝒟,x]

Plot[{PDF[𝒟,x],CDF[𝒟,x],SurvivalFunction[𝒟,x]},{x,-5,5},PlotPoints->100,Filling->Axis]

Quantile[𝒟, p]

Plot[%,{p,0,1}]

Moment related functions:

Table[Moment[𝒟,n],{n,5}]

Probabilities and Expectations:

NProbability[x^2<x,x\[Distributed]𝒟]

NExpectation[Exp[-Abs[x]],x\[Distributed]𝒟]

Random variates:

RandomVariate[𝒟,5]

Histogram@RandomVariate[𝒟,10^4]

Basics 2D

Work in any number of dimensions:

𝒟=HistogramDistribution[RandomVariate[NormalDistribution[],{10^3,2}]];

All properties work:

Table[Plot3D[df[𝒟,{x,y}],{x,-3,3},{y,-3,3},Exclusions->None],{df,{PDF,CDF,SurvivalFunction}}]

Probabilities & Expectation:

NProbability[x^2<1&&x>-1/2,{x,y}\[Distributed]𝒟]

Random variates:

Histogram3D[RandomVariate[𝒟,10^4]]

Marginal distributions:

𝒟x=MarginalDistribution[𝒟,1]

Plot[PDF[𝒟x,x],{x,-3,3},PlotPoints->100,Filling->Axis]

Nonparametric Distributions – SmoothKernelDistribution__

(* 
𝒟ref=ComponentMixtureDistribution[{3/4,1/4},{NormalDistribution[],NormalDistribution[3/2,1/3]}];
datap=ListLinePlot[RandomNumber[𝒟ref,500],Ticks->None,ImageSize->Tiny];
𝒟data=SmoothKernelDistribution[RandomNumber[𝒟ref,10^5]];
distp=Plot[{0,PDF[𝒟data,x]},{x,-3,3},Filling->Axis,Ticks->None,ImageSize->Tiny];
*)


FramedStatement[HoldForm[SmoothKernelDistribution[,…]⟹ ]]

(* 𝒟ref=ComponentMixtureDistribution[{1/3,2/3},{BinormalDistribution[{-1,-1},{1,1},1/2],BinormalDistribution[{1,1},{1,1},-1/2]}];
datap=Graphics[{PointSize[Tiny],Point@RandomNumber[𝒟ref,10^3]},Ticks->None,ImageSize->Tiny];
𝒟dist=SmoothKernelDistribution[RandomNumber[𝒟ref,10^5]];
distp=Plot3D[Evaluate@PDF[𝒟dist,{x,y}],{x,-3,3},{y,-3,3},PlotStyle->Lighter[ColorData[1,2],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,2],0.8],Lighting->"Neutral",Mesh->None,PlotPoints->50,Evaluated->True,ImageSize->Tiny,Ticks->None]; *)

FramedStatement[HoldForm[SmoothKernelDistribution[,…]⟹ ]]

Basics 1D

Suppose we have data coming from the following underlying distribution:

𝒟ref=MixtureDistribution[{3/4,1/4},{NormalDistribution[],NormalDistribution[3/2,1/3]}];

Plot[PDF[𝒟ref,x],{x,-3,3}]

data=RandomVariate[𝒟ref,(10^6)];

Smooth kernel estimate (∑_(i=1)^n)K[((x-(Subscript[x, i]))/h)] where K is the smoothing kernel and h is the bandwidth:

𝒟data=SmoothKernelDistribution[data];

Plot[{PDF[𝒟ref,x],PDF[𝒟data,x]},{x,-3,3}]

Resulting distribution functions:

Table[Plot[f[𝒟data,x],{x,-3,3},PlotLabel->f],{f,{PDF,CDF,SurvivalFunction,HazardFunction}}]

Plot[Quantile[𝒟data,p],{p,0,1}]

Moment related properties:

Table[Moment[𝒟data,n],{n,4}]

Probabilities & Expectations:

NProbability[(x^2)<1,x\[Distributed]𝒟data]

NExpectation[Abs[x-Mean[𝒟data]],x\[Distributed]𝒟data]

Random variates:

RandomVariate[𝒟data,5]

Histogram[RandomVariate[𝒟data,10^5]]

Applications 1D

Compute the temperature distribution for your location:

data=WeatherData[loc=FindGeoLocation[],"Temperature",{{2001,1,1},{2010,1,1}},"NonMetricValue"];
data=Extract[data,Position[data,x_/;NumericQ[x]]]//N;

𝒟=SmoothKernelDistribution[data];

Show[Plot[PDF[𝒟,x],{x,-10,110},Frame->True,Filling->Axis,PlotLabel->Row[{"Temperature Distribution: Lat ",loc[[1]]," Lon ",loc[[2]]}],FrameLabel->{"°F","Density"},Axes->False,PlotRange->All,ColorFunction->Function[{x,y},ColorData["TemperatureMap"][x]]],Plot[PDF[𝒟,x],{x,-10,110},PlotStyle->Black]]

Basics 2D

Work in any number of dimensions:

𝒟ref=MixtureDistribution[{1/3,2/3},{BinormalDistribution[{-1,-1},{1,1},1/2],BinormalDistribution[{1,1},{1,1},-1/2]}];

Plot3D[Evaluate@PDF[𝒟ref,{x,y}],{x,-3,3},{y,-3,3},PlotPoints->35]

data=RandomVariate[𝒟ref,(10^4)];

Get a smooth density distribution:

𝒟data=SmoothKernelDistribution[data];

Plot3D[Evaluate@PDF[𝒟data,{x,y}],{x,-3,3},{y,-3,3},PlotPoints->35]

Moment related:

Table[Moment[𝒟data,{i,j}],{i,3},{j,3}]

Probabilities & Expectations:

NProbability[x<0&&y>2,{x,y}\[Distributed]𝒟data]

Random variates:

RandomVariate[𝒟data,5]

Applications 2D

Old Faithful geyser data: {waiting time, duration}

𝒟faithful=SmoothKernelDistribution[ExampleData[{"Statistics","OldFaithful"}]];

Plot3D[Evaluate@PDF[𝒟faithful,{x,y}],{x,1,5.5},{y,35,105},PlotRange->All,ColorFunction->"Pastel"]

Find the probability that you are waiting more then 3 minutes the duration is less then 80 seconds:

NProbability[waiting>3&&duration<80,{waiting,duration}\[Distributed]𝒟faithful]//Quiet

Formula Distributions

(*
𝒟=ProbabilityDistribution[(Piecewise[{{(E^(x λ)) p λ, x<0}, {(E^(-x λ)) (1-p) λ, x≥0}, {0, True}}]),{x,-∞,∞},Assumptions->λ>0∧0<p<1];
Block[{p=1/3},Plot[Evaluate@Table[PDF[𝒟,x],{λ,{0.5,1}}],{x,-3,3},Filling->Axis,Ticks->None,ImageSize->Tiny]];
*)


FramedStatement[HoldForm[ProbabilityDistribution[(Piecewise[{{(E^(x λ)) p λ, x<0}, {(E^(-x λ)) (1-p) λ, x≥0}}]),…]⟹ ]]

Basics

Define a new distribution by giving the PDF:

𝒟=ProbabilityDistribution[(((Sqrt[2])/π))(1/(1+(x^4))),{x,-∞,∞}]

Behaves just like any other distribution:

Moment[𝒟,2]

PDF[𝒟,x]

Probability[X^2<3,X\[Distributed]𝒟]

Can include parameters, in this case for a two-sided exponential distribution:

𝒟=ProbabilityDistribution[(Piecewise[{{(E^(x λ)) p λ, x<0}, {(E^(-x λ)) (1-p) λ, x≥0}}]),{x,-∞,∞},Assumptions->λ>0∧0<p<1];

CDF[𝒟,x]

Probability[-1/2<X<1/2,X\[Distributed]𝒟]

Applications

Define a reflected gamma distribution:

ReflectedGammaDistribution[α_]=ProbabilityDistribution[Piecewise[{{PDF[GammaDistribution[α,1],x]/2,x≥0},{PDF[GammaDistribution[α,1],-x]/2,x<0}}],{x,-∞,∞},Assumptions->α>0];

pdf=PDF[ReflectedGammaDistribution[α],x]

cdf=CDF[ReflectedGammaDistribution[α],x]

CellDisplay@Manipulate[Plot[{(Piecewise[{{(1/2) ((Piecewise[{{(((E^(-x)) (x^(-1+α)))/(Gamma[α])), x>0}, {0, True}}])), x≥0}, {(1/2) ((Piecewise[{{(((E^x) (((-x))^(-1+α)))/(Gamma[α])), -x>0}, {0, True}}])), x<0}, {0, True}}]),(Piecewise[{{(1/2), x==0}, {((Gamma[α,-x])/(2 Gamma[α])), x<0}, {((2 Gamma[α]-Gamma[α,x])/(2 Gamma[α])), True}}])},{x,-7,7},PlotRange->{0,1},ImageSize->Large],{{α,1.5},0.02,3}]

Define a Gumbel bivariate exponential distribution:

𝒟=ProbabilityDistribution[Exp[-2(x+y)](Exp[x+y]+α(Exp[x]-2)(Exp[y]-2)),{x,0,∞},{y,0,∞},Assumptions->-1<α<1];

CDF[𝒟,{x,y}]

Probability[(X^2)Y<1,{X,Y}\[Distributed]𝒟]

N[%]

Manipulate[Plot3D[(E^(-2 (x+y))) ((E^(x+y))+(-2+(E^x)) (-2+(E^y)) α),{x,0,2},{y,0,2},MeshFunctions->{#3&},MeshShading->{Green,None,Yellow,None},PlotPoints->35,ImageSize->Large,PlotRange->{0,1}],{α,-0.98,0.98}]

Derived Distributions

FramedLargeStatement@Grid[{{TransformedDistribution,"functions of random variables"},
{OrderDistribution,"(n^th) smallest element of random vectors"},
{CensoredDistribution,"interval censored distribution"},
{TruncatedDistribution,"interval truncated distribution"},
{MixtureDistribution,"finite mixture of distributions"},
{ParameterMixtureDistribution,"distribution mixture of distributions"},
{ProductDistribution,"independ multivariate components"},
{CopulaDistribution,"multivariate distribution from kernel and marginals"},
{MarginalDistribution,"lower dimensional distribution"}},
Dividers->{None,{None,None,True,None,True,None,True}},
Alignment->{Left,Baseline}]

Derived Distributions

OrderDistribution

(*
dist=Plot[{0,PDF[NormalDistribution[],x]},{x,-3,3},Filling->Axis,Ticks->None,ImageSize->Tiny];
odistp=Plot[{0,PDF[OrderDistribution[{NormalDistribution[],2},2],x]},{x,-3,3},Filling->Axis,Ticks->None,ImageSize->Tiny];
*)


FramedStatement[HoldForm[OrderDistribution[{,2},2]⟹ ]]

Basics

Lets define a generalization of Min and Max called Nth:

RankedMin[{3,1,2,4},2]

Lets look at the data from picking out the 2nd of 5 elements:

data=Map[RankedMin[#,2]&,RandomVariate[LaplaceDistribution[0,1],{10^3,5}]];

hist=Histogram[data,Automatic,"PDF"]

We can compute the distribution exactly:

PDF[OrderDistribution[{LaplaceDistribution[0,1],5},2],x]

Show[hist,Plot[%,{x,-4,1}]]

These are the other ones:

Plot[Evaluate@Table[PDF[OrderDistribution[{LaplaceDistribution[0,1],5},i],x],{i,5}],{x,-2,2},Filling->Axis]

In some cases OrderDistribution will automaticaly simplify:

OrderDistribution[{UniformDistribution[],n},k]

OrderDistribution[{ExponentialDistribution[λ],n},1]

OrderDistribution[{GumbelDistribution[α,β],n},1]

TransformedDistribution

(*
distp=Plot[{0,CDF[UniformDistribution[{0,1}],x]},{x,0,2},Filling->Axis,Ticks->None,ImageSize->Tiny];
𝒟1=TransformedDistribution[(X^2),X∈UniformDistribution[{0,1}]];
tdistp=Plot[Evaluate@{0,CDF[𝒟1,x]},{x,0,2},Filling->Axis,Ticks->None,ImageSize->Tiny];
*)


FramedStatement[HoldForm[TransformedDistribution[g[x],x\[Distributed]]⟹ ]]

Basics

TransformedDistribution – Distribution of f[x] where x\[Distributed]𝒟

𝒟=TransformedDistribution[(x^2),x\[Distributed]UniformDistribution[{0,1}]]

Get all the properties:

PDF[𝒟,x]//FullSimplify

{Plot[Evaluate@PDF[𝒟,x],{x,0,2},Filling->Axis],Plot[Evaluate@CDF[𝒟,x],{x,0,2},Filling->Axis,PlotRange->{0,1}]}

Random variates:

Histogram@RandomVariate[𝒟,(10^4)]

Many special relations for distributions:

TransformedDistribution[x+y,x\[Distributed]PoissonDistribution[(Subscript[μ, 1])]∧y\[Distributed]PoissonDistribution[(Subscript[μ, 2])]]

TransformedDistribution[(x^2)+(y^2),x\[Distributed]NormalDistribution[]∧y\[Distributed]NormalDistribution[]]

TransformedDistribution[Exp[x],x\[Distributed]NormalDistribution[μ,σ]]

There are surprising simplifications that may happen in several steps:

TransformedDistribution[(x^2),{x,y}\[Distributed]DirichletDistribution[{1,2,3}]]

Derived Distributions

TruncatedDistribution

(*
distp=Plot[{0,PDF[NormalDistribution[0,1],x]},{x,-3,3},Filling->Axis,Ticks->None,PlotRange->{0,0.5},ImageSize->Tiny];
truncp=Plot[{0,PDF[TruncatedDistribution[{-1,2},NormalDistribution[0,1]],x]},{x,-3,3},Filling->Axis,Ticks->None,PlotRange->{0,0.5},ImageSize->Tiny];
*)



FramedStatement[HoldForm[TruncatedDistribution[{(Subscript[x, min]),(Subscript[x, max])},]⟹ ]]

Basics

Truncate the values to an interval — or condition on the interval:

𝒟=TruncatedDistribution[{-1,1},NormalDistribution[0,1]]

PDF[𝒟,x]

Plot[{PDF[𝒟,x],PDF[NormalDistribution[0,1],x]},{x,-2,2},Filling->Axis]

CensoredDistribution

(*
distp=Plot[{0,CDF[NormalDistribution[0,1],x]},{x,-3,3},Filling->Axis,Ticks->None,ImageSize->Tiny];
censp=Plot[Evaluate@{0,CDF[CensoredDistribution[{-1,∞},NormalDistribution[0,1]],x]},{x,-3,3},Exclusions->None,Filling->Axis,Ticks->None,ImageSize->Tiny];
*)

FramedStatement[HoldForm[CensoredDistribution[{(Subscript[x, min]),(Subscript[x, max])},]⟹]]

Basics

A left-censored distribution:

𝒟=CensoredDistribution[{a,∞},NormalDistribution[0,1]]

CDF[𝒟,x]

Manipulate[Plot[%,{x,-3,3},Exclusions->None,Filling->Axis],{{a,-1},-3,3}]

A right-censored distribution:

𝒟=CensoredDistribution[{-∞,b},NormalDistribution[0,1]]

CDF[𝒟,x]

Manipulate[Plot[%,{x,-3,3},Exclusions->None,Filling->Axis],{{b,1},-3,3}]

MixtureDistribution

(*
mix1p=Plot[{0,PDF[NormalDistribution[1,1/2],x]},{x,-1,3},Filling->Axis,Ticks->None,ImageSize->Tiny,PlotRange->{0,2.5}];
mix2p=Plot[{0,PDF[NormalDistribution[2,1/6],x]},{x,-1,3},Filling->Axis,Ticks->None,ImageSize->Tiny,PlotRange->{0,2.5}];
mixp=Plot[{0,PDF[MixtureDistribution[{3/5,2/5}, {NormalDistribution[1,1/2],NormalDistribution[2,1/6]}],x]},{x,-1,3},Filling->Axis,Ticks->None,ImageSize->Tiny,PlotRange->{0,2.5}];
*)



FramedStatement[HoldForm[MixtureDistribution[{(Subscript[λ, 1]),(Subscript[λ, 2])},{,}]⟹ ]]

Basics

Define a new distribution as a component mixture of two normal distributions:

𝒟=MixtureDistribution[{3/5,2/5}, {NormalDistribution[1,1/2],NormalDistribution[2,1/6]}]

PDF[𝒟,x]

Plot[{PDF[𝒟,x],PDF[NormalDistribution[1,1/2],x],PDF[NormalDistribution[2,1/6],x]},{x,-1,3},Filling->Axis,PlotRange->All]

Applications

Can produce rich variety of distribution shapes with relatively few components:

Gaussian = NormalDistribution[];
SkewedUnimodal = MixtureDistribution[{1/5,1/5,3/5},{NormalDistribution[],NormalDistribution[1/2,2/3],NormalDistribution[13/12,5/9]}];
StrongSkew = MixtureDistribution[ConstantArray[1/8,8],Table[NormalDistribution[3((2/3)^i - 1),(2/3)^i],{i,0,7}]];
KurtoticUnimodal = MixtureDistribution[{2/3,1/3},{NormalDistribution[],NormalDistribution[0,1/10]}];
Outlier = MixtureDistribution[{1/10,9/10},{NormalDistribution[],NormalDistribution[0,1/10]}];
Bimodal = MixtureDistribution[{1/2,1/2},{NormalDistribution[-1,2/3],NormalDistribution[1,2/3]}];
SeparatedBimodal = MixtureDistribution[{1/2,1/2},{NormalDistribution[-3/2,1/2],NormalDistribution[3/2,1/2]}];
SkewedBimodal = MixtureDistribution[{3/4,1/4},{NormalDistribution[],NormalDistribution[3/2,1/3]}];
Trimodal = MixtureDistribution[{9/20,9/20,1/10},{NormalDistribution[-6/5,3/5],NormalDistribution[6/5,3/5],NormalDistribution[0,1/4]}];
Claw = MixtureDistribution[{1/2,1/10,1/10,1/10,1/10,1/10},Join[{NormalDistribution[]},Table[NormalDistribution[i/2 - 1, 1/10],{i,0,4}]]];
DoubleClaw = MixtureDistribution[Join[{49/100,49/100},ConstantArray[1/350,7]],Join[{NormalDistribution[-1,2/3],NormalDistribution[1,2/3]},Table[NormalDistribution[(i-3)/2,1/100],{i,0,6}]]];
AsymmetricClaw = MixtureDistribution[Join[{1/2},Table[(2^(1-i))/31,{i,-2,2}]],Join[{NormalDistribution[]},Table[NormalDistribution[i + 1/2, (2^-i)/10],{i,-2,2}]]];
AsymmetricDoubleClaw = MixtureDistribution[{46/100,46/100,1/300,1/300,1/300,7/300,7/300,7/300},Join[Join[Table[NormalDistribution[2i - 1,2/3],{i,0,1}],Table[NormalDistribution[-i/2,1/100],{i,1,3}]],Table[NormalDistribution[i/2,7/100],{i,1,3}]]];
SmoothComb = MixtureDistribution[Table[(2^(5-i))/63,{i,0,5}],Table[NormalDistribution[(65-(96(1/2)^i))/21,(32/63)/2^i],{i,0,5}]];
DiscreteComb = MixtureDistribution[{2/7,2/7,2/7,1/21,1/21,1/21},Join[Table[NormalDistribution[(12i - 15)/7,2/7],{i,0,2}],Table[NormalDistribution[2i/7,1/21],{i,8,10}]]];
dists={Gaussian,SkewedUnimodal,StrongSkew,KurtoticUnimodal,Outlier,Bimodal,SeparatedBimodal,SkewedBimodal,Trimodal,Claw,DoubleClaw,AsymmetricClaw,AsymmetricDoubleClaw,SmoothComb,DiscreteComb};
plots=Table[Plot[PDF[𝒟,x],{x,-3,3},Filling->Axis,FillingStyle->Directive[Opacity[0.4],Lighter[ColorData[1,2],0.5]],PlotStyle->ColorData[1,2],Exclusions->None,Evaluated->True,PlotPoints->100,PlotRange->All,AxesOrigin->{0,0}],{𝒟,dists}];
GraphicsGrid[Partition[plots,3],ImageSize->Large]

Derived Distributions

MarginalDistribution

(*
distp=Show[Plot3D[CDF[BinormalDistribution[1/2],{x,y}],{x,-2,2},{y,-2,2},PlotStyle->Lighter[ColorData[1,2],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,2],0.8],MeshFunctions->{#3&},Lighting->"Neutral",Ticks->None,PlotRange->All,ImageSize->Tiny],ParametricPlot3D[{{u,2,CDF[BinormalDistribution[1/2],{u,2}]},{2,u,CDF[BinormalDistribution[1/2],{2,u}]}},{u,-2,2},PlotStyle->{Yellow,Green}]/.Line[pts_]->Tube[pts,0.05],PlotRange->{0,1}];
margp=Plot[{0,CDF[BinormalDistribution[1/2],{u,2}]},{u,-2,2},Filling->Axis,Ticks->None,ImageSize->Tiny];
*)


FramedStatement[HoldForm[MarginalDistribution[,…]⟹ ]]

Basics

Suppose we have a multivariate distribution built from univariate:

𝒟=ProductDistribution[NormalDistribution[μ,σ],ExponentialDistribution[λ],ChiSquareDistribution[ν]];

These are 1D marginal distributions:

MarginalDistribution[𝒟,1]

These are 2D marginal distributions:

MarginalDistribution[𝒟,{1,2}]

In some case marginals are known to be other special distributions:

MarginalDistribution[BinormalDistribution[ρ],1]

MarginalDistribution[MultinomialDistribution[n, {p,q}], 1]

In other cases there is no special result, but it still a valid distribution:

𝒟=MarginalDistribution[MultinomialDistribution[4, {1/2,1/3,1/6}], {1,2}]

PDF[𝒟, {x,y}]

DiscretePlot3D[%,{x,0,4},{y,0,4},ExtentSize->0.5,PlotRange->All]

The distribution function together with its marginals:

Table[MarginalDistributionPlot3D[DirichletDistribution[{2,3,2}],df,{{0,1},{0,1}}],{df,{PDF,CDF,SurvivalFunction}}]

ProductDistribution

(*
ind1p=Plot[{0,PDF[NormalDistribution[0,1],x]},{x,-2,2},Filling->Axis,Ticks->None,ImageSize->Tiny];
ind2p=Plot[{0,PDF[LaplaceDistribution[0,1],x]},{x,-2,2},Filling->Axis,Ticks->None,ImageSize->Tiny];
indp=Plot3D[{0,PDF[ProductDistribution[NormalDistribution[0,1],LaplaceDistribution[0,1]],{x,y}]},{x,-2,2},{y,-2,2},PlotStyle->Lighter[ColorData[1,2],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,2],0.8],MeshFunctions->{#3&},Lighting->"Neutral",Ticks->None,PlotRange->All,ImageSize->Tiny];
*)



FramedStatement[HoldForm[ProductDistribution[,]⟹ ]]

Basics

ProductDistribution – distribution for independent random variables

𝒟=ProductDistribution[NormalDistribution[0,1],LaplaceDistribution[0,1]]

PDF[𝒟,{x,y}]

Plot3D[%,{x,-2,2},{y,-2,2}]

CopulaDistribution

(*
subkp=Plot3D[Evaluate@CDF[CopulaDistribution["Maximal",{UniformDistribution[{0,1}],UniformDistribution[{0,1}]}],{x,y}],{x,0,1},{y,0,1},PlotStyle->Lighter[ColorData[1,2],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,2],0.8],MeshFunctions->{#3&},Lighting->"Neutral",Ticks->None,ImageSize->Tiny];
marg1p=Plot[{0,CDF[ExponentialDistribution[2],x]},{x,-1,3},Filling->Axis,Ticks->None,ImageSize->Tiny];
marg2p=Plot[{0,CDF[NormalDistribution[],x]},{x,-3,3},Filling->Axis,Ticks->None,ImageSize->Tiny];
copp=Plot3D[Evaluate@CDF[CopulaDistribution["Maximal",{ExponentialDistribution[2],NormalDistribution[]}],{x,y}],{x,-1,3},{y,-3,3},PlotStyle->Lighter[ColorData[1,2],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,2],0.8],MeshFunctions->{#3&},Lighting->"Neutral",Ticks->None,PlotPoints->35,ImageSize->Tiny];
*)



FramedStatement[HoldForm[CopulaDistribution[,{,}]⟹ ]]

Basics

CopulaDistribution – a multivariate CDF C[(Subscript[F, 1]),…,(Subscript[F, n])] with subkernel C and marginal distributions (Subscript[F, i])[(Subscript[x, i])]

Fmax=Refine[CDF[CopulaDistribution["Maximal",{UniformDistribution[{0,1}],UniformDistribution[{0,1}]}],{x,y}],0<x<1&&0<y<1]

Fmin=Refine[CDF[CopulaDistribution["Minimal",{UniformDistribution[{0,1}],UniformDistribution[{0,1}]}],{x,y}],0<x<1&&0<y<1]

minmax=Plot3D[{Fmax,Fmin},{x,0,1},{y,0,1},PlotLabel->"Frechet Lower/Upper Bound",ImageSize->Medium,PlotStyle->Opacity[0.5],Mesh->10,MeshFunctions->{#3&},MeshShading->{Yellow,None,Green,None}]

Lots of different copula families to model dependency structures:

ParallelTable[Plot3D[Evaluate@CDF[CopulaDistribution[c,{UniformDistribution[{0,1}],UniformDistribution[{0,1}]}],{x,y}],{x,0,1},{y,0,1},PlotLabel->c],{c,{"Minimal","Maximal","Product",{"Frank",1/2},{"Clayton",5},{"GumbelHougaard",2},{"FGM",1/3},{"AMH",1/2},{"Binormal",1/2}}}]

Use all possible univariate marginal distributions:

ParallelTable[Plot3D[Evaluate@CDF[CopulaDistribution[{"Frank",1/2},{UniformDistribution[{0,1}],d}],{x,y}],{x,0,1},{y,-3,3},PlotLabel->d],{d,{UniformDistribution[{0,1}],ExponentialDistribution[2],NormalDistribution[0,1],LaplaceDistribution[0,1],GumbelDistribution[1,2]}}]

Statistical Data

Curated Data

Use weather data to get a wind distribution:

data=WeatherData["Chicago", "WindSpeed", {{2006, 1, 1}, {2010, 1,1}, "Day"}];

dist=EstimatedDistribution[data[[All,2]],WeibullDistribution[α,β]]

Show[Histogram[data[[All,2]],Automatic,"PDF"],Plot[PDF[dist,x],{x,0,50},PlotStyle->Thick]]

Efficiency for a GE 1.5MW wind turbine:

turbine=Interpolation[{{0.`,0.`},{0.5`,0.`},{1.`,0.`},{1.5`,0.`},{2.`,0.`},{2.5`,0.`},{3.`,0.`},{3.5`,0.`},{4.`,36.`},{4.5`,66.`},{5.`,104.`},{5.5`,150.`},{6.`,205.`},{6.5`,269.`},{7.`,344.`},{7.5`,428.`},{8.`,528.`},{8.5`,644.`},{9.`,774.`},{9.5`,926.5`},{10.`,1079.`},{10.5`,1211.`},{11.`,1342.`},{11.5`,1401.`},{12.`,1460.`},{12.5`,1477.`},{13.`,1494.`},{13.5`,1500.`},{14.`,1500.`},{14.5`,1500.`},{15.`,1500.`},{15.5`,1500.`},{16.`,1500.`},{16.5`,1500.`},{17.`,1500.`},{17.5`,1500.`},{18.`,1500.`},{18.5`,1500.`},{19.`,1500.`},{19.5`,1500.`},{20.`,1500.`},{20.5`,1500.`},{21.`,1500.`},{21.5`,1500.`},{22.`,1500.`},{22.5`,1500.`},{23.`,1500.`},{23.5`,1500.`},{24.`,1500.`},{24.5`,1500.`},{25.`,1500.`},{25.5`,0.`},{26.`,0.`},{26.5`,0.`},{27.`,0.`},{27.5`,0.`},{28.`,0.`},{28.5`,0.`},{29.`,0.`},{29.5`,0.`},{30.`,0.`}}];

Expected energy production in a year:

NExpectation[24 365 turbine[x],x\[Distributed]dist]

ExampleData

A collection of statistical data sets:

ExampleData["Statistics"]

Length[%]

ExampleData[{"Statistics","AirplaneGlass"},"Properties"]

ExampleData[{"Statistics","AirplaneGlass"},"Description"]

Statistical Visualization

Distribution Shapes

Binned and smooth histograms in 1D:

data1=RandomVariate[NormalDistribution[0,1],200];

{Histogram[data1,Automatic,"PDF"],SmoothHistogram[data1,Filling->Axis]}

Binned and smooth histograms in 2D:

data2=RandomVariate[BinormalDistribution[1/2],1000];

{Histogram3D[data2,Automatic,"PDF"],SmoothHistogram3D[data2]}

{DensityHistogram[data2,Automatic,"PDF"],SmoothDensityHistogram[data2]}

Paired histograms:

data1=RandomVariate[NormalDistribution[],200];
data2=RandomVariate[WeibullDistribution[1,2],200];

PairedHistogram[data1,data2]

Distribution Fits

data1=RandomVariate[UniformDistribution[{0,1}],100];

Quantile-quantile plot:

QuantilePlot[data1,NormalDistribution[]]

Probability-probability plot:

ProbabilityPlot[data1]

Scaled Probability-Probability plot, i.e. Weibull, normal, lognormal plots etc:

ProbabilityScalePlot[data1,"Weibull",GridLines->Automatic]

Distribution Comparisons

Compare the essential quantiles of multiple datasets with box-and-whisker:

data=Table[RandomVariate[NormalDistribution[μ,1],100],{μ,0,4,1}];

BoxWhiskerChart[data]

Or compare the full distributions of multiple datasets:

DistributionChart[data]

Initialization

Display Boxes

CellDisplay: print cell without Out[nn]

CellDisplay[expr_]:=CellPrint[ExpressionCell[expr,"Output",Background->None,ShowCellLabel->False,CellFrame->False]]

FrameDisplay: frame object with standard settings for this document

FrameDisplay[expr_]:=Framed[expr,RoundingRadius->10,Background->Lighter[Gray,0.9],FrameMargins->10,FrameStyle->GrayLevel[0.5]]

Use this for definitions:

DefinitionGrid[{{duv_,dmv_},{cuv_,cmv_}}]:=
FrameDisplay@Grid[
{{SpanFromLeft,Text[Style["Univariate",FontFamily->"Helvetica"]],Text[Style["Multivariate",FontFamily->"Helvetica"]]},
{Text[Style["Discrete",FontFamily->"Helvetica"]],TraditionalForm@duv,TraditionalForm@dmv},
{Text[Style["Continuous",FontFamily->"Helvetica"]],TraditionalForm@cuv,TraditionalForm@cmv}},
Spacings->{2,4},
Dividers->{{False,True,True,False},{False,True,True,False}}
]

FormattedDefinitionGrid[{{duv_,dmv_},{cuv_,cmv_}}]:=
FrameDisplay@Style[Grid[
{{SpanFromLeft,Text[Style["Univariate",FontFamily->"Helvetica"]],Text[Style["Multivariate",FontFamily->"Helvetica"]]},
{Text[Style["Discrete",FontFamily->"Helvetica"]],duv,dmv},
{Text[Style["Continuous",FontFamily->"Helvetica"]],cuv,cmv}},
Spacings->{2,4},
Dividers->{{False,True,True,False},{False,True,True,False}}
],"Text"]

For moment definitions:

MomentDefinitionGrid[{{uvm_,uvgf_},{mvm_,mvgf_}}]:=
FrameDisplay@
Style[Grid[
{{SpanFromLeft,Text[Style["Moment",FontFamily->"Helvetica"]],Text[Style["Generating Function",FontFamily->"Helvetica"]]},
{Text[Style["Univariate",FontFamily->"Helvetica"]],TraditionalForm@uvm,TraditionalForm@uvgf},
{Text[Style["Multivariate",FontFamily->"Helvetica"]],TraditionalForm@mvm,TraditionalForm@mvgf}},
Spacings->{2,4},
Dividers->{{False,True,True,False},{False,True,True,False}}
],"Text"]

For moments:

MomentGrid[list_,type_]:=
Grid[Transpose@{list},Alignment->{Left,Automatic},Background->{None,{{Lighter[Blend[Switch[type,Moment,{Blue,Green},CentralMoment,{Red,Orange},FactorialMoment,{Magenta,Cyan},Cumulant,{Purple,Red}]],.8],White}}},Dividers->{None,{Darker[Gray,.6],{False},Darker[Gray,.6]}},ItemSize->30]

MomentTabView[dist_,n_Integer?Positive]:=
CellDisplay@TabView[{"Moment"->MomentGrid[Table[Moment[dist,k],{k,n}],Moment],"CentralMoment"->MomentGrid[Table[CentralMoment[dist,k],{k,n}],CentralMoment],"FactorialMoment"->MomentGrid[Table[FactorialMoment[dist,k],{k,n}],FactorialMoment],"Cumulant"->MomentGrid[Table[Cumulant[dist,k],{k,n}],Cumulant]}]

For moment generating functions:

MGFTabView[dist_]:=
CellDisplay@TabView[{"MGF"->MomentGeneratingFunction[dist,t],"CMGF"->CentralMomentGeneratingFunction[dist,t],"FMGF"-> FactorialMomentGeneratingFunction[dist,t],"CGF"->CumulantGeneratingFunction[dist,t],"CF"->CharacteristicFunction[dist,t]}]

FramedStatement:

FramedStatement[expr_]:=
CellDisplay@FrameDisplay@Text[Style[expr,Large,FontFamily->"Helvetica"]]

FramedLargeStatement[expr_]:=
CellDisplay@FrameDisplay@Text[Style[expr,FontSize->24,FontFamily->"Helvetica"]]

FramedMediumStatement[expr_]:=
CellDisplay@FrameDisplay@Text[Style[expr,FontSize->16,FontFamily->"Helvetica"]]

Used for enumerating properties:

FormulaGallery[forms_List] := CellDisplay@TraditionalForm@Grid[forms, Dividers->All,Spacings->{1,2}, Alignment ->{Center,Center}, BaseStyle->{FontFamily->"Verdana"},Background->{None,{{Hue[.6,.15,.9],GrayLevel[.9]}}},FrameStyle->Directive[Thick,White]]

DistributionPlot, DistributionPlot3D

Notes/Todos: 
• Should probably store default plot range for each distribution. The quantile based heuristic works fine in 90% of cases only. 
• Could consider adding support for multiple distributions. 
• Could consider using a PlotLayout option which is either “CommonXAxis” or “Grid” ... 
• Could consider controls for what elements to plot.

DataPlot1D: Useful ways of displaying point clouds.

ListPlot1D[data_List]:=DataPlot1D[data,"JitteredPoints"]

ListPlot1D[data_List,"Points",opts___?OptionQ]:=
Graphics[{PointSize[Tiny],Point@Developer`ToPackedArray@Transpose[{data,ConstantArray[0.,Length[data]]}]},FilterOptions[Graphics,opts],AspectRatio->1/10]

ListPlot1D[data_List,"JitteredPoints",opts___?OptionQ]:=
Graphics[{PointSize[Tiny],Point@Transpose[{data,RandomReal[{-1,1},Length[data]]}]},FilterOptions[Graphics,opts],AspectRatio->1/10]

ListPlot1D[data_List,"Lines",opts___?OptionQ]:=
Graphics[{Thickness[Tiny],Line@Developer`ToPackedArray@Table[{{data[[i]],-1.},{data[[i]],1}},{i,Length[data]}]},FilterOptions[Graphics,opts],AspectRatio->1/10]

DistribubtionPlot related functions:

ClearAll[DistributionPlot,DistributionPlot3D]

Statistics`Library`ContinuousUnivariateDistributionQ[ExponentialDistribution[2]];

Utilities:

BivariateContinuousDistributionQ[dist_]:=Statistics`Library`ContinuousMultivariateDistributionQ[dist]&&Statistics`Library`DistributionDimensionality[dist]==2

BivariateDiscreteDistributionQ[dist_]:=Statistics`Library`DiscreteMultivariateDistributionQ[dist]&&Statistics`Library`DistributionDimensionality[dist]==2

FilterOptions:

FilterOptions[s_Symbol,opts___?OptionQ]:=Sequence@@FilterRules[{opts},Options[s]]

DistributeDefinitions[FilterOptions]

Defaults:

DistributionPlot[dist_,opts___?OptionQ]/;Statistics`Library`UnivariateDistributionQ[dist]:=DistributionPlot[dist,Quantile[dist,{0.05,0.95}],opts]

Continuous univariate distributions:

DistributionPlot[dist_,{xmin_,xmax_},opts___?OptionQ]/;Statistics`Library`ContinuousUnivariateDistributionQ[dist]:=
Module[{x,pdf,cdf,sf,hf},
{pdf,cdf,sf,hf}= 
{
Plot[PDF[dist,x],{x,xmin,xmax},PlotRange->{0,Automatic},PlotLabel->"PDF",PlotStyle->ColorData[1,1],Filling->Axis,FillingStyle->Lighter[ColorData[1,1],0.7],Evaluate@FilterOptions[Plot,opts]],Plot[CDF[dist,x],{x,xmin,xmax},PlotRange->{0,1},PlotLabel->"CDF",PlotStyle->ColorData[1,2],Filling->Axis,FillingStyle->Lighter[ColorData[1,2],0.7],Evaluate@FilterOptions[Plot,opts]],Plot[SurvivalFunction[dist,x],{x,xmin,xmax},PlotRange->{0,1},PlotLabel->"SF",PlotStyle->ColorData[1,3],Filling->Axis,FillingStyle->Lighter[ColorData[1,3],0.7],Evaluate@FilterOptions[Plot,opts]],
Plot[HazardFunction[dist,x],{x,xmin,xmax},PlotRange->{0,Automatic},PlotLabel->"HF",PlotStyle->ColorData[1,4],Filling->Axis,FillingStyle->Lighter[ColorData[1,4],0.7],Evaluate@FilterOptions[Plot,opts]]
};

Framed[GraphicsGrid[{{pdf,cdf},{sf,hf}},ImageSize->Large],RoundingRadius->10,Background->Lighter[Gray,0.9]]
]

Discrete univariate distributions:

DistributionPlot[dist_,{xmin_,xmax_},opts___?OptionQ]/;Statistics`Library`DiscreteUnivariateDistributionQ[dist]:=
Block[{k,x,kmin=Floor[xmin],kmax=Ceiling[xmax],pdf,cdf,sf,hf},
{pdf,cdf,sf,hf}=
{
DiscretePlot[PDF[dist,k],{k,kmin,kmax},PlotRange->{{xmin,xmax},Automatic},PlotLabel->"PDF",PlotStyle->ColorData[1,1],Filling->Axis,FillingStyle->Lighter[ColorData[1,1],0.7],ExtentSize->0.5,Evaluate@FilterOptions[DiscretePlot,opts]],
DiscretePlot[CDF[dist,k],{k,kmin,kmax},PlotRange->{0,1},PlotLabel->"CDF",PlotStyle->ColorData[1,2],Filling->Axis,FillingStyle->Lighter[ColorData[1,2],0.7],ExtentSize->Right,ExtentMarkers->{"Filled","Empty"},Evaluate@FilterOptions[Plot,opts]],
DiscretePlot[SurvivalFunction[dist,k],{k,kmin,kmax},PlotRange->{0,1},PlotLabel->"SF",PlotStyle->ColorData[1,3],Filling->Axis,FillingStyle->Lighter[ColorData[1,3],0.7],ExtentSize->Right,ExtentMarkers->{"Filled","Empty"},Evaluate@FilterOptions[Plot,opts]],
DiscretePlot[HazardFunction[dist,k],{k,kmin,kmax},PlotRange->{{xmin,xmax},Automatic},PlotLabel->"HF",PlotStyle->ColorData[1,4],Filling->Axis,FillingStyle->Lighter[ColorData[1,4],0.7],ExtentSize->0.5,Evaluate@FilterOptions[DiscretePlot,opts]]
};

Framed[GraphicsGrid[{{pdf,cdf},{sf,hf}},ImageSize->Large],RoundingRadius->10,Background->Lighter[Gray,0.9]]
]

Continuous bivariate distributions:

DistributionPlot3D[dist_,{xmin_,xmax_},{ymin_,ymax_},opts___?OptionQ]/;BivariateContinuousDistributionQ[dist]:=
Module[{x,y,pdf,cdf,sf,hf},
{pdf,cdf,sf,hf}=
{Plot3D[PDF[dist,{x,y}],{x,xmin,xmax},{y,ymin,ymax},PlotRange->{0,Automatic},PlotLabel->"PDF",PlotStyle->Lighter[ColorData[1,1],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,1],0.8],MeshFunctions->{#3&},Lighting->"Neutral",Evaluate@FilterOptions[Plot3D,opts]],
Plot3D[CDF[dist,{x,y}],{x,xmin,xmax},{y,ymin,ymax},PlotRange->{0,1},PlotLabel->"CDF",PlotStyle->Lighter[ColorData[1,2],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,2],0.8],MeshFunctions->{#3&},Lighting->"Neutral",Evaluate@FilterOptions[Plot3D,opts]],
Plot3D[SurvivalFunction[dist,{x,y}],{x,xmin,xmax},{y,ymin,ymax},PlotRange->{0,1},PlotLabel->"SF",PlotStyle->Lighter[ColorData[1,3],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,3],0.8],MeshFunctions->{#3&},Lighting->"Neutral",Evaluate@FilterOptions[Plot3D,opts]],
Plot3D[HazardFunction[dist,{x,y}],{x,xmin,xmax},{y,ymin,ymax},PlotRange->{0,Automatic},PlotLabel->"HF",PlotStyle->Lighter[ColorData[1,4],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,4],0.8],MeshFunctions->{#3&},Lighting->"Neutral",Evaluate@FilterOptions[Plot3D,opts]]};

Framed[GraphicsGrid[{{pdf,cdf},{sf,hf}},ImageSize->Large],RoundingRadius->10,Background->Lighter[Gray,0.9]]
]

Discrete bivariate distributions:

DistributionPlot3D[dist_,{xmin_,xmax_},{ymin_,ymax_},opts___?OptionQ]/;BivariateDiscreteDistributionQ[dist]:=
Block[{x,y,n,m,nmin=Floor[xmin],nmax=Ceiling[xmax],mmin=Floor[ymin],mmax=Ceiling[ymax],pdf,cdf,sf,hf},
{pdf,cdf,sf,hf}=
{DiscretePlot3D[PDF[dist,{n,m}],{n,nmin,nmax},{m,mmin,mmax},PlotLabel->"PDF",PlotStyle->Lighter[ColorData[1,1],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,1],0.8],Lighting->"Neutral",ExtentSize->0.5,Evaluate@FilterOptions[DiscretePlot3D,opts]],
DiscretePlot3D[CDF[dist,{x,y}],{x,xmin,xmax},{y,ymin,ymax},PlotRange->{0,1},PlotLabel->"CDF",PlotStyle->Lighter[ColorData[1,2],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,2],0.8],ExtentSize->Right,Lighting->"Neutral",Evaluate@FilterOptions[Plot3D,opts]],
DiscretePlot3D[SurvivalFunction[dist,{x,y}],{x,xmin,xmax},{y,ymin,ymax},PlotRange->{0,1},PlotLabel->"SF",PlotStyle->Lighter[ColorData[1,3],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,3],0.8],ExtentSize->Right,Lighting->"Neutral",Evaluate@FilterOptions[Plot3D,opts]],
DiscretePlot3D[HazardFunction[dist,{n,m}],{n,nmin,nmax},{m,mmin,mmax},PlotLabel->"HF",ExtentSize->0.5,PlotStyle->Lighter[ColorData[1,4],0.5],Filling->Axis,FillingStyle->Lighter[ColorData[1,4],0.8],Lighting->"Neutral",Evaluate@FilterOptions[DiscretePlot3D,opts]]};

Framed[GraphicsGrid[{{pdf,cdf},{sf,hf}},ImageSize->Large],RoundingRadius->10,Background->Lighter[Gray,0.9]]
]

Other Visualization

A particular style of comparing histograms with density functions:

CompareHistogram3D[dist_,{xmin_,xmax_,dx_},{ymin_,ymax_,dy_},df_:PDF,opts___?OptionQ]:=
Module[{hist,dplot},
hist=Histogram3D[RandomVariate[dist,10^4],{{xmin,xmax,dx},{ymin,ymax,dy}},ToString[df],PlotRange->{{xmin,xmax},{ymin,ymax},All},ChartBaseStyle->Opacity[0.5]];
dplot=Plot3D[df[dist,{x,y}],{x,xmin,xmax},{y,ymin,ymax},Mesh->{Range[xmin,xmax,dx],Range[ymin,ymax,dy]},PlotRange->All,MeshStyle->Gray,PlotStyle->Hue[.15,.7,.8]];
Show[hist,dplot]
];

Showing the marginals on walls:

MarginalDistributionPlot3D[dist_,df_:PDF,{{xmin_,xmax_},{ymin_,ymax_}},opts___?OptionQ]:=
Module[{gxy,gx,gy},
gxy=Plot3D[Evaluate@df[dist,{x,y}],{x,xmin,xmax},{y,ymin,ymax},Mesh->None,PlotRange->All];
gx=Plot[Evaluate@df[MarginalDistribution[dist,1],x],{x,xmin,xmax},Filling->Axis,FillingStyle->Red];
gy=Plot[Evaluate@df[MarginalDistribution[dist,2],y],{y,ymin,ymax},Filling->Axis,FillingStyle->Green];
Graphics3D[{First[gxy],First[gx]/.GraphicsComplex[pts_,rest__]:>GraphicsComplex[pts/.{x_,y_}->{x,ymin,y},rest],First[gy]/.GraphicsComplex[pts_,rest__]:>GraphicsComplex[pts/.{x_,y_}->{xmin,x,y},rest]},BoxRatios->{1,1,0.4}]
]

Distribution Converter

PDF[ExtremeValueDistribution[α,β],x]== PDF[MaxStableDistribution[α,β,0],x]//Simplify

Simplify[PDF[FrechetDistribution[α,β],x]== PDF[MaxStableDistribution[β,β/α,1/α],x],DistributionParameterAssumptions[FrechetDistribution[α,β]]]

Simplify[PDF[WeibullDistribution[α,β],x]==PDF[MinStableDistribution[β,β/α,-1/α],x],DistributionParameterAssumptions[WeibullDistribution[α,β]]]

PDF[GumbelDistribution[α,β],x]==PDF[MinStableDistribution[α,β,0],x]//Simplify

DistributionConvert:

DistributionConvert[ExtremeValueDistribution[α_,β_]]:=MaxStableDistribution[α,β,0]

DistributionConvert[FrechetDistribution[α_,β_]]:=MaxStableDistribution[β,β/α,1/α]

DistributionConvert[WeibullDistribution[α_,β_]]:=MinStableDistribution[β,β/α,-1/α]

DistributionConvert[GumbelDistribution[μ_,β_]]:=MinStableDistribution[α,β,0]

Extra:

//FramedStatement

Wrapup

SetDirectory[NotebookDirectory[]];

$Line=0;