--- title: "EstimatedDistribution" language: "en" type: "Symbol" summary: "EstimatedDistribution[data, dist] estimates the parametric distribution dist from data. EstimatedDistribution[data, dist, {{p, p0}, {q, q0}, ...}] estimates the parameters p, q, ... with starting values p0, q0, .... EstimatedDistribution[data, dist, idist] estimates distribution dist with starting values taken from the instantiated distribution idist." keywords: - distribution fitting - distribution estimation - parameter estimation - maximum likelihood - ml - mle - ml estimate - ml estimator - maximum likelihood estimator - maximum likelihood estimate - profile maximum likelihood - profile maximum likelihood estimate - method of moments - mom - method of moments estimate - method of moments estimator canonical_url: "https://reference.wolfram.com/language/ref/EstimatedDistribution.html" source: "Wolfram Language Documentation" related_guides: - title: "Probability & Statistics with Quantities" link: "https://reference.wolfram.com/language/guide/ProbabilityWithQuantities.en.md" - title: "Reliability" link: "https://reference.wolfram.com/language/guide/Reliability.en.md" - title: "Probability & Statistics" link: "https://reference.wolfram.com/language/guide/ProbabilityAndStatistics.en.md" - title: "Random Variables" link: "https://reference.wolfram.com/language/guide/RandomVariables.en.md" - title: "Survival Analysis" link: "https://reference.wolfram.com/language/guide/SurvivalAnalysis.en.md" - title: "Statistical Data Analysis" link: "https://reference.wolfram.com/language/guide/Statistics.en.md" - title: "Optimization" link: "https://reference.wolfram.com/language/guide/Optimization.en.md" - title: "Scientific Data Analysis" link: "https://reference.wolfram.com/language/guide/ScientificDataAnalysis.en.md" - title: "Robust Descriptive Statistics" link: "https://reference.wolfram.com/language/guide/RobustDescriptiveStatistics.en.md" - title: "Descriptive Statistics" link: "https://reference.wolfram.com/language/guide/DescriptiveStatistics.en.md" - title: "Tabular Modeling" link: "https://reference.wolfram.com/language/guide/TabularModeling.en.md" - title: "Time Series Processing" link: "https://reference.wolfram.com/language/guide/TimeSeries.en.md" - title: "Actuarial Computation" link: "https://reference.wolfram.com/language/guide/ActuarialComputation.en.md" related_functions: - title: "FindDistributionParameters" link: "https://reference.wolfram.com/language/ref/FindDistributionParameters.en.md" - title: "FindDistribution" link: "https://reference.wolfram.com/language/ref/FindDistribution.en.md" - title: "LogLikelihood" link: "https://reference.wolfram.com/language/ref/LogLikelihood.en.md" - title: "Moment" link: "https://reference.wolfram.com/language/ref/Moment.en.md" - title: "DistributionFitTest" link: "https://reference.wolfram.com/language/ref/DistributionFitTest.en.md" - title: "EmpiricalDistribution" link: "https://reference.wolfram.com/language/ref/EmpiricalDistribution.en.md" - title: "SmoothKernelDistribution" link: "https://reference.wolfram.com/language/ref/SmoothKernelDistribution.en.md" - title: "QuantityDistribution" link: "https://reference.wolfram.com/language/ref/QuantityDistribution.en.md" --- # EstimatedDistribution EstimatedDistribution[data, dist] estimates the parametric distribution dist from data. EstimatedDistribution[data, dist, {{p, p0}, {q, q0}, …}] estimates the parameters p, q, … with starting values p0, q0, …. EstimatedDistribution[data, dist, idist] estimates distribution dist with starting values taken from the instantiated distribution idist. ## Details and Options * ``EstimatedDistribution`` returns the distribution ``dist`` with parameter estimates inserted for any non-numeric values. * The ``data`` must be a list of possible outcomes from the given distribution ``dist``. * The distribution ``dist`` can be any parametric univariate, multivariate, or derived distribution with unknown parameters. * The following options can be given: | | | | | ------------------- | ------------------- | ------------------------------------------- | | AccuracyGoal | Automatic | the accuracy sought | | ParameterEstimator | "MaximumLikelihood" | what parameter estimator to use | | PrecisionGoal | Automatic | the precision sought | | WorkingPrecision | Automatic | the precision used in internal computations | * The following basic settings can be used for ``ParameterEstimator``: | | | | -------------------------- | ------------------------------------ | | "MaximumLikelihood" | maximize the log‐likelihood function | | "MethodOfMoments" | match raw moments | | "MethodOfCentralMoments" | match central moments | | "MethodOfCumulants" | match cumulants | | "MethodOfFactorialMoments" | match factorial moments | * The maximum likelihood method attempts to maximize the log-likelihood function $\sum _i \log \left(f\left(x_i;\theta \right)\right)$, where $\theta$ are the distribution parameters and $f(x;\theta )$ is the PDF of the distribution. * The method of moments solves $\hat{\mu }_i=\mu _i(\theta )$, $\hat{\mu }_j=\mu _j(\theta )$, $\ldots$, where $\hat{\mu }_i$ is the $i$$$^{\text{th}}$$ sample moment and $\mu _i(\theta )$ is the $i$$$^{\text{th}}$$ moment of the distribution, with parameters $\theta$. * Method-of-moment-based estimators may not satisfy all restrictions on parameters. --- ## Examples (47) ### Basic Examples (3) Obtain the maximum likelihood parameter estimates, assuming a gamma distribution: ```wl In[1]:= data = RandomVariate[\[ScriptCapitalD] = GammaDistribution[1, 2], 10 ^ 3]; In[2]:= e\[ScriptCapitalD] = EstimatedDistribution[data, GammaDistribution[α, β]] Out[2]= GammaDistribution[0.929364, 2.03594] ``` Visually compare the PDFs for the original and estimated distributions: ```wl In[3]:= Plot[{PDF[\[ScriptCapitalD], x], PDF[e\[ScriptCapitalD], x]}, {x, 0, 10}, PlotLegends -> {"\[ScriptCapitalD]", "e\[ScriptCapitalD]"}] Out[3]= [image] ``` Obtain the method of moments estimates: ```wl In[4]:= EstimatedDistribution[data, GammaDistribution[α, β], ParameterEstimator -> "MethodOfMoments"] Out[4]= GammaDistribution[0.872599, 2.16839] ``` --- Estimate parameters for a multivariate distribution: ```wl In[1]:= data = RandomVariate[MultivariatePoissonDistribution[5, {4, 3, 10}], 20]; In[2]:= EstimatedDistribution[data, MultivariatePoissonDistribution[μ, {Subscript[μ, 1], Subscript[μ, 2], Subscript[μ, 3]}]] Out[2]= MultivariatePoissonDistribution[4.33788, {5.96212, 4.71212, 11.4621}] ``` --- Estimated parameters from data with quantities: ```wl In[1]:= data = RandomVariate[NormalDistribution[Quantity[200, "Centimeters"], Quantity[7, "Centimeters"]], 20] Out[1]= StructuredArray[QuantityArray, {20}, StructuredArray`StructuredData[QuantityArray, {203.8978123846134, 197.4365777389725, 200.40858430012597, 193.6052040986552, 193.28819094692227, 193.61173625013265, 205.90084715548764, 199.76596469092476, ... .69271988807867, 203.06606250223138, 207.21369511284902, 192.29925813725887, 196.71839444871966, 194.6865134520807, 202.17995515095512, 207.23336004021073, 193.45828547669208, 203.13072769409234, 186.33164621969425}, "Centimeters", {{1}}]] In[2]:= EstimatedDistribution[data, NormalDistribution[m, s]] Out[2]= QuantityDistribution[NormalDistribution[199.317, 6.04233], "Centimeters"] ``` ### Scope (15) #### Basic Uses (5) Estimate both parameters for a binomial distribution: ```wl In[1]:= bdata = {23, 16, 24, 17, 16, 17, 18, 17, 19, 23, 19, 13, 19, 18, 22, 21, 19, 22, 19, 23}; In[2]:= EstimatedDistribution[bdata, BinomialDistribution[n, p]] Out[2]= BinomialDistribution[33, 0.583333] ``` Estimate ``p``, assuming ``n`` is known: ```wl In[3]:= EstimatedDistribution[bdata, BinomialDistribution[50, p]] Out[3]= BinomialDistribution[50, 0.385] ``` Estimate ``n``, assuming ``p`` is known: ```wl In[4]:= EstimatedDistribution[bdata, BinomialDistribution[n, 2 / 5]] Out[4]= BinomialDistribution[48, (2/5)] ``` --- Get the distribution with maximum likelihood parameter estimate for a particular family: ```wl In[1]:= vals = RandomReal[ChiSquareDistribution[20], 1000]; In[2]:= dist = EstimatedDistribution[vals, ChiSquareDistribution[n]] Out[2]= ChiSquareDistribution[19.9805] ``` Check goodness of fit by comparing a histogram of the data and the estimate's PDF: ```wl In[3]:= Show[Histogram[vals, Automatic, "ProbabilityDensity"], Plot[PDF[dist, x], {x, 0, 50}]] Out[3]= [image] ``` Perform goodness-of-fit tests with null distribution ``dist`` : ```wl In[4]:= DistributionFitTest[vals, dist, {"TestDataTable", All}] ``` Out[4]= Statistic P-Value Anderson-Darling 1.07836 0.318592 Cramér-von Mises 0.173818 0.324504 Kolmogorov-Smirnov 0.0308104 0.292558 Kuiper 0.047481 0.145308 Pearson \[Chi]^2 39.616 0.137931 Watson U^2 0.138382 0.130331 Perform tests correcting for estimation of the parameter: ```wl In[5]:= DistributionFitTest[vals, ChiSquareDistribution[n], {"TestDataTable", All}] ``` Out[5]= Statistic P-Value Anderson-Darling 1.07836 0.0507442 Cramér-von Mises 0.173818 0.0417188 Kolmogorov-Smirnov 0.0308104 0.0417647 Kuiper 0.047481 0.0968893 Pearson \[Chi]^2 39.616 0.112518 Watson U^2 0.138382 0.0787759 --- Estimate parameters by maximizing the log‐likelihood: ```wl In[1]:= data = RandomVariate[GeometricDistribution[.4], 100]; In[2]:= dist = EstimatedDistribution[data, GeometricDistribution[p]] Out[2]= GeometricDistribution[0.37594] In[3]:= phat = dist[[1]] Out[3]= 0.37594 ``` Plot the log‐likelihood function to visually check that the solution is optimal: ```wl In[4]:= Show[Plot[LogLikelihood[GeometricDistribution[p], data], {p, 0, 1}], Graphics[{PointSize[Large], Red, Point[{phat, LogLikelihood[GeometricDistribution[phat], data]}]}]] Out[4]= [image] ``` --- Visualize a log‐likelihood surface to find rough values for the parameters: ```wl In[1]:= data = BlockRandom[SeedRandom[100];RandomVariate[WeibullDistribution[4, 5], 50]]; In[2]:= ContourPlot[LogLikelihood[WeibullDistribution[α, β], data], {α, 2, 6}, {β, 3, 8}] Out[2]= [image] ``` Supply those rough values as starting values for the estimation: ```wl In[3]:= EstimatedDistribution[data, WeibullDistribution[α, β], WeibullDistribution[3, 5.5]] Out[3]= WeibullDistribution[3.63331, 5.42278] ``` --- Estimate the normal approximation of Poisson data: ```wl In[1]:= data = RandomVariate[PoissonDistribution[20], 100]; In[2]:= EstimatedDistribution[data, NormalDistribution[μ, σ]] Out[2]= NormalDistribution[20.11, 4.62579] ``` Obtain the estimate to 20 digits: ```wl In[3]:= EstimatedDistribution[data, NormalDistribution[μ, σ], WorkingPrecision -> 20] Out[3]= NormalDistribution[20.110000000000000000, 4.6257864196263968884] ``` #### Univariate Parametric Distributions (2) Estimate parameters for a continuous distribution: ```wl In[1]:= data = {0.8, 0.98, 0.66, 0.09, 0.03, 0.16, 0.88, 0.32, 0.56, 0.3, 0.82, 0.34, 0.93, 0.33, 0.25, 0.74, 0.3, 0.74, 0.16, 0.02}; In[2]:= EstimatedDistribution[data, BetaDistribution[α, β]] Out[2]= BetaDistribution[0.832619, 0.921159] ``` Compare empirical and distribution quantiles: ```wl In[3]:= QuantilePlot[data, %] Out[3]= [image] ``` --- Estimate parameters for a discrete distribution: ```wl In[1]:= data = {16, 1, 2, 5, 3, 5, 1, 6, 4, 14, 2, 9, 1, 1, 2, 2, 8, 1, 3, 1}; In[2]:= EstimatedDistribution[data, LogSeriesDistribution[θ]] Out[2]= LogSeriesDistribution[0.914521] ``` #### Multivariate Parametric Distributions (2) Estimate parameters for a discrete multivariate distribution: ```wl In[1]:= mdata = RandomVariate[MultinomialDistribution[20, {1 / 3, 1 / 6, 1 / 10, 2 / 5}], 1000]; In[2]:= EstimatedDistribution[mdata, MultinomialDistribution[n, {Subscript[p, 1], Subscript[p, 2], Subscript[p, 3], Subscript[p, 4]}]] Out[2]= MultinomialDistribution[20, {0.33375, 0.1655, 0.0964, 0.40435}] ``` --- Estimate parameters for a continuous multivariate distribution: ```wl In[1]:= bndata = RandomVariate[BinormalDistribution[{1, 2}, {1 / 3, 4}, 3 / 4], 1000]; In[2]:= dist = EstimatedDistribution[bndata, BinormalDistribution[{Subscript[μ, 1], Subscript[μ, 2]}, {Subscript[σ, 1], Subscript[σ, 2]}, ρ]] Out[2]= BinormalDistribution[{0.984464, 1.80515}, {0.334629, 3.82977}, 0.727583] ``` Compare the difference between the original and estimated PDFs: ```wl In[3]:= Plot3D[PDF[BinormalDistribution[{1, 2}, {1 / 3, 4}, 3 / 4], {x, y}] - PDF[dist, {x, y}], {x, 0, 2}, {y, -10, 16}, PlotRange -> All] Out[3]= [image] ``` #### Derived Distributions (6) Estimate parameters for a truncated normal: ```wl In[1]:= \[ScriptCapitalD]o = TruncatedDistribution[{-1, 2}, NormalDistribution[1, 1 / 3]]; In[2]:= SeedRandom[1]; data = RandomVariate[\[ScriptCapitalD]o, 100]; In[3]:= \[ScriptCapitalD]e = EstimatedDistribution[data, TruncatedDistribution[{-1, 2}, NormalDistribution[μ, σ]]] Out[3]= TruncatedDistribution[{-1, 2}, NormalDistribution[0.998497, 0.338917]] ``` Compare original and estimated distribution: ```wl In[4]:= Plot[{PDF[\[ScriptCapitalD]o, x], PDF[\[ScriptCapitalD]e, x]}, {x, -1, 3}, Filling -> Axis, PlotLegends -> {"\[ScriptCapitalD]o", "\[ScriptCapitalD]e"}] Out[4]= [image] ``` --- Estimate parameters for a constructed distribution: ```wl In[1]:= dist = ProbabilityDistribution[λ Exp[λ(1 - x)], {x, 1, Infinity}, Assumptions -> λ > 0] Out[1]= ProbabilityDistribution[E^(1 - \[FormalX]) λ λ, {\[FormalX], 1, ∞}, Assumptions -> λ > 0] In[2]:= EstimatedDistribution[{1.4, 5.1, 1.7, 1.6, 1.1, 3.9, 2.2, 1.3, 2., 1.5}, dist] Out[2]= ProbabilityDistribution[0.847458 E^0.847458 (1 - \[FormalX]), {\[FormalX], 1, ∞}] ``` --- Estimate parameters for a product distribution: ```wl In[1]:= data = RandomVariate[ProductDistribution[BetaDistribution[2, 3], LaplaceDistribution[10, 1]], 100]; In[2]:= EstimatedDistribution[data, ProductDistribution[BetaDistribution[α, β], LaplaceDistribution[μ, σ]]] Out[2]= ProductDistribution[BetaDistribution[1.71741, 3.11346], LaplaceDistribution[9.78866, 1.13259]] ``` --- Estimate parameters for a copula distribution: ```wl In[1]:= \[ScriptCapitalD]o = CopulaDistribution[{"Frank", 1 / 2}, {GammaDistribution[2, 3], ChiSquareDistribution[10]}]; In[2]:= data = RandomVariate[\[ScriptCapitalD]o, 100]; In[3]:= \[ScriptCapitalD]e = EstimatedDistribution[data, CopulaDistribution[{"Frank", c}, {GammaDistribution[α, β], ChiSquareDistribution[n]}]] Out[3]= CopulaDistribution[{"Frank", 0.693803}, {GammaDistribution[1.87505, 3.36944], ChiSquareDistribution[10.252]}] ``` Compare original and estimated CDFs: ```wl In[4]:= Plot3D[CDF[\[ScriptCapitalD]o, {x, y}] - CDF[\[ScriptCapitalD]e, {x, y}], {x, 0, 25}, {y, 0, 35}, PlotRange -> All] Out[4]= [image] ``` --- Estimate parameters for a component mixture: ```wl In[1]:= data = RandomVariate[MixtureDistribution[{1 / 3, 2 / 3}, {GammaDistribution[2, 3], NormalDistribution[1, 1 / 4]}], 100]; In[2]:= EstimatedDistribution[data, MixtureDistribution[{p, 1 - p}, {GammaDistribution[α, β], NormalDistribution[μ, σ]}]]//Quiet Out[2]= MixtureDistribution[{0.278526, 0.721474}, {GammaDistribution[2.1646, 2.4859], NormalDistribution[0.981615, 0.287674]}] ``` Estimate the mixture probabilities assuming the component distributions are known: ```wl In[3]:= EstimatedDistribution[data, MixtureDistribution[{p, 1 - p}, {GammaDistribution[2, 3], NormalDistribution[1, 1 / 4]}]]//Quiet Out[3]= MixtureDistribution[{0.286951, 0.713049}, {GammaDistribution[2, 3], NormalDistribution[1, (1/4)]}] ``` --- Estimate parameters for quantity distribution in specified units: ```wl In[1]:= data = RandomVariate[NormalDistribution[Quantity[200, "Centimeters"], Quantity[7, "Centimeters"]], 100] Out[1]= StructuredArray[QuantityArray, {100}, StructuredArray`StructuredData[QuantityArray, CompressedData["«1152»"], "Centimeters", {{1}}]] In[2]:= EstimatedDistribution[data, QuantityDistribution[NormalDistribution[m, s], "Meters"]] Out[2]= QuantityDistribution[NormalDistribution[1.99687, 0.0600811], "Meters"] ``` ### Options (4) #### ParameterEstimator (3) Estimate parameters by matching central moments: ```wl In[1]:= data = RandomVariate[ParetoDistribution[2, 8, 5], 1000]; In[2]:= EstimatedDistribution[data, ParetoDistribution[α, β, μ], ParameterEstimator -> "MethodOfCentralMoments"] Out[2]= ParetoDistribution[5.83513, 21.1659, 4.9963] ``` Other moment‐based methods typically give similar results: ```wl In[3]:= EstimatedDistribution[data, ParetoDistribution[α, β, μ], ParameterEstimator -> "MethodOfCentralMoments"] Out[3]= ParetoDistribution[5.83513, 21.1659, 4.9963] In[4]:= EstimatedDistribution[data, ParetoDistribution[α, β, μ], ParameterEstimator -> "MethodOfFactorialMoments"] Out[4]= ParetoDistribution[4.35826, 16.3229, 5.00123] ``` --- Estimate parameters based on default moments: ```wl In[1]:= data = RandomVariate[NormalDistribution[2, 3], 1000]; In[2]:= EstimatedDistribution[data, NormalDistribution[μ, σ], ParameterEstimator -> "MethodOfMoments"] Out[2]= NormalDistribution[1.92707, 3.05945] ``` Estimate parameters from the first and fourth moments: ```wl In[3]:= EstimatedDistribution[data, NormalDistribution[μ, σ], ParameterEstimator -> {"MethodOfMoments", "MomentOrders" -> {1, 4}}] Out[3]= NormalDistribution[1.92707, 3.09387] ``` --- Obtain the maximum likelihood estimates using the default method: ```wl In[1]:= data = RandomVariate[BetaDistribution[2, 3], 1000]; In[2]:= EstimatedDistribution[data, BetaDistribution[α, β]] Out[2]= BetaDistribution[1.9438, 2.77066] ``` Use ``FindMaximum`` to obtain the estimates: ```wl In[3]:= EstimatedDistribution[data, BetaDistribution[α, β], ParameterEstimator -> {"MaximumLikelihood", Method -> "FindMaximum"}] Out[3]= BetaDistribution[1.9438, 2.77066] ``` Use ``EvaluationMonitor`` to extract the points sampled: ```wl In[4]:= {params, {points}} = Reap[EstimatedDistribution[data, BetaDistribution[α, β], ParameterEstimator -> {"MaximumLikelihood", Method -> {"FindMaximum", EvaluationMonitor :> Sow[{α, β}]}}]]; ``` Visualize the sequences of sampled $\alpha$ and $\beta$ values: ```wl In[5]:= ListLogPlot[Transpose[points], PlotLegends -> {"α", "β"}] Out[5]= [image] ``` #### WorkingPrecision (1) Use machine precision for continuous parameters by default: ```wl In[1]:= data = RandomVariate[GumbelDistribution[2, 5], 100, WorkingPrecision -> 50]; In[2]:= EstimatedDistribution[data, GumbelDistribution[α, β]] Out[2]= GumbelDistribution[1.77748, 4.78591] ``` Obtain a higher-precision result: ```wl In[3]:= EstimatedDistribution[data, GumbelDistribution[α, β], WorkingPrecision -> 25] Out[3]= GumbelDistribution[1.777476803556768655217093, 4.785906855228969096863159] ``` ### Applications (14) #### Estimation of Similarly Shaped Distributions (1) Model lognormal distributed data with a gamma distribution: ```wl In[1]:= lnorm = LogNormalDistribution[2, 0.3]; In[2]:= sample = RandomVariate[lnorm, 10 ^ 4]; In[3]:= edist = EstimatedDistribution[sample, GammaDistribution[α, β]] Out[3]= GammaDistribution[11.2557, 0.686048] In[4]:= Show[Histogram[sample, 20, "ProbabilityDensity"], Plot[PDF[edist, x], {x, 0, 20}, PlotStyle -> Thick]] Out[4]= [image] ``` Compare the distributions of the simulation and estimated distributions: ```wl In[5]:= Plot[{PDF[lnorm, x], PDF[edist, x]}, {x, 0, 25}] Out[5]= [image] ``` #### Accident Claims (1) The number of accident claims per policy per year from an insurance company: ```wl In[1]:= accidentsPerPolicyCounts = {81714, 11306, 1618, 250, 40, 7};; In[2]:= accidentsPerPolicy = WeightedData[Range[0, 5], accidentsPerPolicyCounts] Out[2]= WeightedData[Automatic, {{0, 1, 2, 3, 4, 5}, {81714, 11306, 1618, 250, 40, 7}}] ``` Model the data by a logarithmic series distribution since most policies have at most one claim: ```wl In[3]:= edist = EstimatedDistribution[WeightedData[Range[0, 5] + 1, accidentsPerPolicyCounts], LogSeriesDistribution[θ]] Out[3]= LogSeriesDistribution[0.255459] In[4]:= Show[Histogram[accidentsPerPolicy, {-0.5, 5.5, 1}, "PDF"], DiscretePlot[PDF[edist, x + 1], {x, 0, 5}, PlotStyle -> PointSize[Medium]]] Out[4]= [image] ``` #### Word Lengths in Different Languages (1) Get word length data for several languages: ```wl In[1]:= languages = {"Arabic", "English", "Finnish", "French", "Hebrew", "Hindi", "Italian", "Russian", "Spanish"}; In[2]:= worddata = Table[StringLength /@ DictionaryLookup[{l, All}], {l, languages}]; ``` Model the word lengths for each language as binomially distributed: ```wl In[3]:= binom = Table[EstimatedDistribution[i, BinomialDistribution[n, p], ParameterEstimator -> {"MaximumLikelihood", Method -> {"FindRoot", MaxIterations -> 1000}}], {i, worddata}] Out[3]= {BinomialDistribution[20, 0.312457], BinomialDistribution[38, 0.220887], BinomialDistribution[91, 0.139018], BinomialDistribution[30, 0.320026], BinomialDistribution[17, 0.39756], BinomialDistribution[33, 0.181917], BinomialDistribution[26, 0.370831], BinomialDistribution[31, 0.261445], BinomialDistribution[29, 0.300194]} ``` Compare the actual and estimated distributions: ```wl In[4]:= Partition[Table[Show[Histogram[worddata[[i]], {Range[25] - 1 / 2}, "ProbabilityDensity", PlotLabel -> languages[[i]]], DiscretePlot[PDF[binom[[i]], x], {x, 0, 25}, PlotRange -> All, PlotStyle -> PointSize[.025]]], {i, Length[languages]}], 3]//Grid Out[4]= | | | | | ------- | ------- | ------- | | [image] | [image] | [image] | | [image] | [image] | [image] | | [image] | [image] | [image] | ``` #### Text Frequency (1) The word count in a text follows a Zipf distribution: ```wl In[1]:= text = ExampleData[{"Text", "OriginOfSpecies"}, "Words"]; In[2]:= wordCount = Tally[text][[All, 2]]; ``` Fit a ``ZipfDistribution`` to the word frequency data: ```wl In[3]:= edist = EstimatedDistribution[wordCount, ZipfDistribution[ρ]] Out[3]= ZipfDistribution[0.759346] ``` Compare the frequency histogram with the estimated distribution: ```wl In[4]:= Show[Histogram[wordCount, {0.5, 20.5, 1}, "ProbabilityDensity"], DiscretePlot[PDF[edist, x], {x, 0, 20}, PlotStyle -> PointSize[Medium]]] Out[4]= [image] ``` #### Earthquake Magnitudes (1) ``EstimatedDistribution`` can be used with constructs like ``MixtureDistribution`` to create multimodal models: ```wl In[1]:= magnitudes = Select[ExampleData[{"Statistics", "USEarthquakes"}], #[[1]] ≥ 1935&][[All, 7]]; ``` The magnitudes of earthquakes in the United States in the selected years have two modes: ```wl In[2]:= h = Histogram[magnitudes, 20, "ProbabilityDensity"] Out[2]= [image] ``` Fit distribution from possible mixtures of one ``NormalDistribution`` with another: ```wl In[3]:= edist = EstimatedDistribution[magnitudes, MixtureDistribution[{p, 1 - p}, {NormalDistribution[a, b], NormalDistribution[c, d]}]]//Quiet Out[3]= MixtureDistribution[{0.174389, 0.825611}, {NormalDistribution[4.52858, 0.0503323], NormalDistribution[5.35465, 0.808865]}] ``` Compare the histogram to the PDF of the estimated distribution: ```wl In[4]:= Show[h, Plot[PDF[edist, x], {x, 0, 10}, PlotStyle -> Thick, PlotRange -> All]] Out[4]= [image] ``` Find the probability of an earthquake of magnitude 7 or higher: ```wl In[5]:= Probability[x ≥ 7, x\[Distributed]edist] Out[5]= 0.0173118 ``` Find the mean earthquake magnitude: ```wl In[6]:= Mean[edist] Out[6]= 5.21059 ``` Simulate magnitudes of the next 30 earthquakes: ```wl In[7]:= ListPlot[RandomVariate[edist, 30], Filling -> Axis] Out[7]= [image] ``` #### Wind Speed Analysis (1) Model monthly maximum wind speeds in Boston: ```wl In[1]:= maxWinds = WeatherData["Boston", "MaxWindSpeed", {{1950, 1, 1}, {2009, 12, 31}, "Month"}, "Value"]//QuantityMagnitude; ``` Fit the data to a ``RayleighDistribution`` : ```wl In[2]:= edist1 = EstimatedDistribution[maxWinds, RayleighDistribution[θ]] Out[2]= RayleighDistribution[37.907] ``` An ``ExtremeValueDistribution`` : ```wl In[3]:= edist2 = EstimatedDistribution[maxWinds, ExtremeValueDistribution[a, b]] Out[3]= ExtremeValueDistribution[47.2118, 9.02774] ``` Compare the empirical quantiles and those for the fitted distributions to see where the models deviate from the data: ```wl In[4]:= Row[Table[QuantilePlot[maxWinds, i, PlotLabel -> Head[i]], {i, {edist1, edist2}}]] Out[4]= [image][image] ``` #### Distribution of Incomes (1) Model incomes at a large state university: ```wl In[1]:= ExampleData[{"Statistics", "UniversitySalaries"}, "ColumnDescriptions"] Out[1]= {"Department", "Estimated percentage of full time employment.", "Projected 2010 salary.", "Campus location"} In[2]:= universityA = Cases[ExampleData[{"Statistics", "UniversitySalaries"}], {dept_, perc_, salary_, campus : "A"} :> {perc, salary}]; In[3]:= salaries = Cases[universityA, {p_ ? Positive, s_ ? Positive} :> s / p]; ``` Assume the salaries are Dagum distributed: ```wl In[4]:= edist = EstimatedDistribution[salaries, DagumDistribution[p, a, b]] Out[4]= DagumDistribution[5.75551, 2.66493, 23588.5] ``` Assume they follow a more general Pareto distribution: ```wl In[5]:= edist2 = EstimatedDistribution[salaries, ParetoDistribution[k, α, γ, μ], ParameterEstimator -> {"MaximumLikelihood", "Method" -> "SimulatedAnnealing"}] Out[5]= ParetoDistribution[36126.8, 0.612621, 0.264385, 7200.] ``` Compare the subtle differences in the estimated distributions: ```wl In[6]:= Show[Histogram[salaries, {0, 200000, 10000}, "PDF"], Plot[{PDF[edist, x], PDF[edist2, x]}, {x, 0, 200000}, PlotStyle -> Thick, PlotLegends -> {DagumDistribution, ParetoDistribution}]] Out[6]= [image] ``` #### Automobile Fuel Efficiency (1) The average city and highway mileage for midsize cars follows a binormal distribution: ```wl In[1]:= midsizeMPG2009 = {{16, 22}, {18, 26}, {17, 25}, {16, 23}, {16, 23}, {14, 19}, {13, 19}, {10, 14}, {10, 17}, {18, 28}, {18, 27}, {17, 25}, {17, 25}, {17, 26}, {17, 26}, {16, 25}, {17, 25}, {15, 22}, {15, 23}, {11, 17}, {11, 17}, {17, 28}, {16, 24}, {16, 25}, {16, 25}, {17, 26}, {18, 26}, {17, 26}, {17, 25}, {17, 26}, {13, 19}, {15, 24}, {17, 26}, {15, 22}, {22, 33}, {22, 30}, {18, 29}, {17, 26}, {26, 34}, {21, 30}, {13, 20}, {19, 27}, {16, 27}, {21, 30}, {13, 20}, {19, 27}, {16, 27}, {24, 30}, {23, 27}, {23, 29}, {21, 25}, {19, 27}, {10, 15}, {9, 16}, {17, 25}, {20, 29}, {20, 28}, {18, 26}, {24, 33}, {25, 33}, {16, 25}, {15, 23}, {22, 32}, {22, 32}, {20, 28}, {23, 30}, {24, 32}, {19, 27}, {19, 26}, {18, 25}, {17, 24}, {16, 24}, {16, 23}, {16, 24}, {16, 23}, {20, 22}, {17, 25}, {20, 29}, {20, 28}, {18, 26}, {17, 24}, {18, 28}, {20, 29}, {21, 30}, {17, 25}, {18, 25}, {23, 32}, {17, 24}, {16, 22}, {15, 22}, {13, 19}, {13, 20}, {20, 27}, {16, 25}, {23, 32}, {23, 31}, {18, 27}, {19, 26}, {35, 33}, {19, 26}, {24, 30}, {21, 29}, {26, 31}, {27, 33}, {24, 32}, {11, 18}, {24, 30}, {24, 32}, {22, 33}, {17, 26}, {26, 34}, {21, 31}, {21, 31}, {19, 28}, {33, 34}, {48, 45}, {19, 29}, {15, 23}, {15, 22}, {16, 25}}; ``` Assume city and highway miles per gallon are normally distributed and correlated: ```wl In[2]:= edist = EstimatedDistribution[midsizeMPG2009, BinormalDistribution[{μ1, μ2}, {σ1, σ2}, ρ]] Out[2]= BinormalDistribution[{18.6167, 26.25}, {5.0186, 4.66949}, 0.905544] ``` Show the distribution of city and highway mileage: ```wl In[3]:= Plot[Evaluate[PDF[MarginalDistribution[edist, #], x]& /@ {1, 2}], {x, 5, 40}, Filling -> Axis, PlotLegends -> {"city", "highway"}] Out[3]= [image] ``` Visualize the joint density with contours on a logarithmic scale: ```wl In[4]:= ContourPlot[Log[PDF[edist, {x, y}]], {x, 0, 50}, {y, 0, 50}, Contours -> 50] Out[4]= [image] ``` #### Earthquake Waiting Times (1) The data contains waiting times in days between serious (magnitude at least 7.5 or over 1000 fatalities) earthquakes worldwide, recorded from 12/16/1902 to 3/4/1977: ```wl In[1]:= earthquakesWaitingTimes = ExampleData[{"Statistics", "EarthquakeWaitingTimes"}]; ``` Model waiting times by an ``ExponentialDistribution`` : ```wl In[2]:= edist = EstimatedDistribution[earthquakesWaitingTimes, ExponentialDistribution[λ]] Out[2]= ExponentialDistribution[0.00228723] ``` Estimate the average and median number of days between major earthquakes: ```wl In[3]:= {Mean[edist], Median[edist]} Out[3]= {437.21, 303.051} ``` #### Earthquake Frequency (1) The number of earthquakes per year can be modeled by ``SinghMaddalaDistribution`` : ```wl In[1]:= ExampleData[{"Statistics", "USEarthquakes"}, "ColumnDescriptions"] Out[1]= {"The year.", "The month.", "The day.", "The UTC time hhmmss.mm.", "The latitude.", "The longitude.", "The magnitude."} In[2]:= earthquakes = Tally[Select[ExampleData[{"Statistics", "USEarthquakes"}], #[[1]] ≥ 1935&][[All, 1]]][[All, 2]] Out[2]= {29, 27, 25, 38, 31, 56, 42, 36, 22, 22, 19, 35, 33, 30, 33, 38, 29, 98, 22, 89, 64, 42, 130, 60, 70, 36, 49, 52, 47, 102, 136, 70, 41, 46, 57, 39, 54, 38, 45, 46, 78, 39, 37, 35, 55, 141, 40, 38, 64, 40, 40, 94, 49, 55, 58} ``` Fit the distribution to the data: ```wl In[3]:= edist = EstimatedDistribution[earthquakes, SinghMaddalaDistribution[q, a, b]] Out[3]= SinghMaddalaDistribution[0.470391, 5.63701, 35.0974] ``` Compare the data histogram with the PDF of the estimated distribution: ```wl In[4]:= Show[Histogram[earthquakes, 15, "ProbabilityDensity"], Plot[PDF[edist, x], {x, 0, 140}, PlotStyle -> Thick]] Out[4]= [image] ``` Find the probability of at least 60 earthquakes in the US in a year: ```wl In[5]:= NProbability[x ≥ 60, x\[Distributed]edist] Out[5]= 0.235937 ``` #### Time between Geyser Eruptions (1) Mixtures can be used to model multimodal data: ```wl In[1]:= ExampleData[{"Statistics", "OldFaithful"}, "ColumnDescriptions"] Out[1]= {"Eruption time in minutes", "Waiting time to next eruption in minutes"} In[2]:= waiting = ExampleData[{"Statistics", "OldFaithful"}][[All, 2]]; ``` A histogram of waiting times for eruptions of the Old Faithful geyser exhibits two modes: ```wl In[3]:= h = Histogram[waiting, 20, "ProbabilityDensity"] Out[3]= [image] ``` Fit a ``MixtureDistribution`` to the data: ```wl In[4]:= mdist = MixtureDistribution[{1 / 3, 2 / 3}, {GammaDistribution[a, b], GammaDistribution[c, d]}]; In[5]:= edist = Quiet@EstimatedDistribution[waiting, mdist, {{a, 85}, {b, .5}, {c, 195}, {d, .4}}] Out[5]= MixtureDistribution[{(1/3), (2/3)}, {GammaDistribution[83.3088, 0.657739], GammaDistribution[193.172, 0.415161]}] ``` Compare the histogram to the PDF of the estimated distribution: ```wl In[6]:= Show[h, Plot[PDF[edist, x], {x, 0, 100}, PlotStyle -> Thick]] Out[6]= [image] ``` Find the probability that the waiting time is over 80 minutes: ```wl In[7]:= NProbability[x > 80, x\[Distributed]edist] Out[7]= 0.0000329977 ``` Simulate waiting times for the next 60 eruptions: ```wl In[8]:= times = RandomVariate[edist, 60]; In[9]:= ListPlot[{times, {{1, m1}, {60, m1}}, {{1, m2}, {60, m2}}}, Joined -> {False, True, True}, Filling -> {1 -> Axis}, AxesOrigin -> {0, 40}] Out[9]= [image] ``` #### Stock Price Distribution (1) Lognormal distribution can be used to model stock prices: ```wl In[1]:= stocks = FinancialData["SBUX", {{2010, 1, 1}, {2014, 1, 1}, "Day"}, "Value"]; ``` Fit the distribution to the data: ```wl In[2]:= edist = EstimatedDistribution[stocks, LogNormalDistribution[μ, σ]] Out[2]= QuantityDistribution[LogNormalDistribution[3.06249, 0.366685], "USDollars"] ``` Observe that the quantiles for the data and distribution match well except for the largest values: ```wl In[3]:= QuantilePlot[stocks, QuantityMagnitude[edist]] Out[3]= [image] ``` #### Water Flow Rates (1) Consider the annual minimum daily flows given in cubic meters per second for the Mahanadi river: ```wl In[1]:= minFlow = ExampleData[{"Statistics", "MahanadiRiverFlow"}]; ``` Model the annual minimum mean daily flows as a ``MinStableDistribution`` : ```wl In[2]:= edist = EstimatedDistribution[minFlow, MinStableDistribution[a, b, c]] Out[2]= MinStableDistribution[2.53019, 1.17611, -0.672677] ``` Compare the histogram of the data to the PDF of the estimated distribution: ```wl In[3]:= Show[Histogram[minFlow, 10, "ProbabilityDensity"], Plot[PDF[edist, x], {x, 0, 7}, PlotStyle -> Thick]] Out[3]= [image] ``` Simulate annual minimum mean daily flows for the next 30 years: ```wl In[4]:= ListLinePlot[RandomVariate[edist, 30]] Out[4]= [image] ``` #### Population Sizes (1) Use a Pareto distribution to model Australian city population sizes: ```wl In[1]:= cities = CityData[{All, "Australia"}]; In[2]:= populations = QuantityMagnitude[CityData[#, "Population"]& /@ cities]; In[3]:= edist = EstimatedDistribution[populations, ParetoDistribution[k, α, γ, μ]] Out[3]= ParetoDistribution[8424.06, 1.11083, 0.920555, 4.] ``` Estimate the probability that a city has a population of at least 10,000 people: ```wl In[4]:= Probability[x ≥ 10000, x\[Distributed]edist] Out[4]= 0.415618 ``` Compute the probability based on the original data: ```wl In[5]:= Probability[x ≥ 10000, x\[Distributed]populations]//N Out[5]= 0.410319 ``` ### Properties & Relations (8) ``EstimatedDistribution`` gives a distribution with parameter estimates inserted: ```wl In[1]:= EstimatedDistribution[{5, 8, 3, 4, 9, 6}, PoissonDistribution[μ]] Out[1]= PoissonDistribution[5.83333] ``` ``FindDistributionParameters`` gives parameter estimates as replacement rules: ```wl In[2]:= FindDistributionParameters[{5, 8, 3, 4, 9, 6}, PoissonDistribution[μ]] Out[2]= {μ -> 5.83333} ``` --- ``EstimatedProcess`` estimates a parametric process: ```wl In[1]:= data = RandomFunction[BernoulliProcess[0.4], {0, 10 ^ 5}]; In[2]:= EstimatedProcess[data, BernoulliProcess[p]] Out[2]= BernoulliProcess[0.397136] ``` ``EstimatedDistribution`` estimates a parametric distribution: ```wl In[3]:= data = RandomVariate[BernoulliDistribution[0.4], 10 ^ 5]; In[4]:= EstimatedDistribution[data, BernoulliDistribution[p]] Out[4]= BernoulliDistribution[0.40179] ``` --- Estimate distribution parameters by maximum likelihood: ```wl In[1]:= data = RandomVariate[ExponentialDistribution[3], 1000]; In[2]:= EstimatedDistribution[data, ExponentialDistribution[λ]] Out[2]= ExponentialDistribution[2.86623] ``` Use ``DistributionFitTest`` to test the quality of the fit: ```wl In[3]:= ℋ = DistributionFitTest[data, ExponentialDistribution[λ], "HypothesisTestData"] Out[3]= HypothesisTestData[DistributionFitTest, {CompressedData["«10658»"], "SampleToNull", 1, Rational[1, 20], {ExponentialDistribution[λ], Automatic}, Automatic, {}}, {"Normality" -> 0, "EqualVariance" -> 0, "Symmetry" -> 0}] ``` Extract the fitted distribution: ```wl In[4]:= ℋ["FittedDistribution"] Out[4]= ExponentialDistribution[2.86623] ``` Obtain a table of relevant test statistics and $p$‐values: ```wl In[5]:= ℋ[{"TestDataTable", All}] ``` Out[5]= Statistic P-Value Anderson-Darling 0.778476 0.218292 Cramér-von Mises 0.115242 0.253422 Kolmogorov-Smirnov 0.0312414 0.0994216 Kuiper 0.0470772 0.10409 Pearson \[Chi]^2 38.272 0.14289 Watson U^2 0.115148 0.139367 --- ``EstimatedDistribution`` estimates parameters in a parametric distribution: ```wl In[1]:= data = RandomVariate[NormalDistribution[0, 1], 50]; In[2]:= edist = EstimatedDistribution[data, NormalDistribution[μ, σ]] Out[2]= NormalDistribution[0.147614, 1.10754] ``` ``SmoothKernelDistribution`` gives a nonparametric kernel density estimate: ```wl In[3]:= skdist = SmoothKernelDistribution[data] Out[3]= DataDistribution["SmoothKernel", {CompressedData["«5546»"], CompressedData["«5610»"], 0.46046436342282254}, 1, 50] ``` Compare the PDFs for the nonparametric and parametric distributions: ```wl In[4]:= Plot[{PDF[skdist, x], PDF[edist, x]}, {x, -3, 3}, PlotStyle -> Thick] Out[4]= [image] ``` Visualize the nonparametric density using ``SmoothHistogram`` : ```wl In[5]:= SmoothHistogram[data, PlotStyle -> Thick] Out[5]= [image] ``` --- ``EstimatedDistribution`` gives a maximum likelihood estimate of parameters: ```wl In[1]:= data = RandomVariate[GammaDistribution[4, 10], 25]; In[2]:= edist = EstimatedDistribution[data, GammaDistribution[α, β]] Out[2]= GammaDistribution[4.60444, 6.85588] ``` Compute the likelihood using ``Likelihood`` : ```wl In[3]:= Likelihood[edist, data] Out[3]= 1.7079520598003813`*^-44 ``` Compute the log‐likelihood using ``LogLikelihood`` : ```wl In[4]:= LogLikelihood[edist, data] Out[4]= -100.778 ``` --- Estimate parameters by matching raw moments: ```wl In[1]:= data = RandomVariate[BetaDistribution[2, 5], 100]; In[2]:= edist = EstimatedDistribution[data, BetaDistribution[α, β], ParameterEstimator -> "MethodOfMoments"] Out[2]= BetaDistribution[2.79602, 6.84447] ``` Compute raw moments from the data using ``Moment`` : ```wl In[3]:= {Moment[data, 1], Moment[data, 2]} Out[3]= {0.290029, 0.103469} ``` Compute the same moments from the estimated distribution: ```wl In[4]:= {Moment[edist, 1], Moment[edist, 2]} Out[4]= {0.290029, 0.103469} ``` --- Estimate parameters for a Weibull distribution: ```wl In[1]:= data = RandomVariate[WeibullDistribution[3, 15], 100]; In[2]:= edist = EstimatedDistribution[data, WeibullDistribution[α, β]] Out[2]= WeibullDistribution[2.59328, 14.4129] ``` Use ``QuantilePlot`` to visualize empirical quantiles versus fitted distribution quantiles: ```wl In[3]:= QuantilePlot[data, edist] Out[3]= [image] ``` Obtain the same visualization when the estimation is done within ``QuantilePlot`` : ```wl In[4]:= QuantilePlot[data, WeibullDistribution[α, β]] Out[4]= [image] ``` --- ``EstimatedDistribution`` ignores time stamps in ``TimeSeries`` and ``EventSeries`` : ```wl In[1]:= ts = TemporalData[TimeSeries, {{{1, 0, 1, 0, 0, 1, 0, 1, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 1, 0, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 1, 1, 1, 0, 0, 1, 0, 0, 0, 1, 1, 1, 0, 1, 1, 1, 1}}, {{0, 99, 1}}, 1, {"Continuous", 1}, {"Discrete", 1}, 1, {ResamplingMethod -> {"Interpolation", InterpolationOrder -> 1}}}, False, 10.1]; In[2]:= EstimatedDistribution[ts, BernoulliDistribution[p]] Out[2]= BernoulliDistribution[0.44] ``` The same as: ```wl In[3]:= EstimatedDistribution[ts["Values"], BernoulliDistribution[p]] Out[3]= BernoulliDistribution[0.44] ``` For ``TemporalData``, all the path structure is ignored: ```wl In[4]:= td = TemporalData[Automatic, {CompressedData["«110034»"], {{0, 1, 0.1}}, 1000, {"Continuous", 1000}, {"Continuous", 1}, 1, {ValueDimensions -> 1, ResamplingMethod -> {"Interpolation", InterpolationOrder -> 1}}}, False, 10.1]; In[5]:= EstimatedDistribution[td, NormalDistribution[a, b]] Out[5]= NormalDistribution[-0.0170957, 0.706739] ``` The same as: ```wl In[6]:= EstimatedDistribution[td["ValueList"]//Flatten, NormalDistribution[a, b]] Out[6]= NormalDistribution[-0.0170957, 0.706739] ``` ### Possible Issues (3) Solutions of method of moment equations can give parameters that are not valid: ```wl In[1]:= EstimatedDistribution[{10, 12, 8, 15}, BinomialDistribution[n, p], ParameterEstimator -> "MethodOfMoments"] Out[1]= BinomialDistribution[27.7397, 0.405556] In[2]:= DistributionParameterQ[%] ``` BinomialDistribution::nnintprm: Parameter 27.73972602739726\` at position 1 in BinomialDistribution[27.7397,0.405556] is expected to be a non-negative integer. ```wl Out[2]= False ``` For a continuous distribution: ```wl In[3]:= data = RandomVariate[NormalDistribution[-1, .1], 10 ^ 3]; In[4]:= EstimatedDistribution[data, NormalDistribution[m, m], ParameterEstimator -> "MethodOfMoments"] Out[4]= NormalDistribution[-0.996756, -0.996756] In[5]:= DistributionParameterQ[%] ``` NormalDistribution::posprm: Parameter -0.996756 at position 2 in NormalDistribution[-0.996756,-0.996756] is expected to be positive. ```wl Out[5]= False ``` --- Good starting values may be needed to obtain a good solution: ```wl In[1]:= data = {4.1, 5.8, 4.3, 5.6, 3.8, 2.5, 3.1, 5.3, 5.4, 3.3}; In[2]:= EstimatedDistribution[data, RiceDistribution[.01, α, β], ParameterEstimator -> "MethodOfMoments"] ``` FindRoot::cvmit: Failed to converge to the requested accuracy or precision within 100 iterations. ```wl Out[2]= RiceDistribution[0.01, 1.847, -3.69954] In[3]:= EstimatedDistribution[data, RiceDistribution[.01, α, β], {{α, 4}, {β, .02}}, ParameterEstimator -> "MethodOfMoments"] Out[3]= RiceDistribution[0.01, 4.4554, 0.108116] ``` --- Good starting values may also result in quicker results: ```wl In[1]:= data = RandomVariate[DirichletDistribution[Range[25]], 10 ^ 4]; In[2]:= EstimatedDistribution[data, DirichletDistribution[Array[g, 25]]]//Timing Out[2]= {1.00126, DirichletDistribution[{0.986669, 2.02751, 2.98826, 3.97894, 5.01635, 5.97562, 7.0092, 8.00742, 8.87519, 9.98504, 11.0175, 12.0049, 12.9686, 13.9941, 15.0089, 16.0247, 16.8965, 17.9206, 18.9795, 20.0015, 21.033, 21.983, 23.0235, 24.0432, 25.0986}]} In[3]:= EstimatedDistribution[data, DirichletDistribution[Array[g, 25]], Transpose[{Array[g, 25], Range[25]}]]//Timing Out[3]= {0.102467, DirichletDistribution[{0.986669, 2.02751, 2.98826, 3.97894, 5.01635, 5.97562, 7.0092, 8.00742, 8.87519, 9.98504, 11.0175, 12.0049, 12.9686, 13.9941, 15.0089, 16.0247, 16.8965, 17.9206, 18.9795, 20.0015, 21.033, 21.983, 23.0235, 24.0432, 25.0986}]} ``` ## See Also * [`FindDistributionParameters`](https://reference.wolfram.com/language/ref/FindDistributionParameters.en.md) * [`FindDistribution`](https://reference.wolfram.com/language/ref/FindDistribution.en.md) * [`LogLikelihood`](https://reference.wolfram.com/language/ref/LogLikelihood.en.md) * [`Moment`](https://reference.wolfram.com/language/ref/Moment.en.md) * [`DistributionFitTest`](https://reference.wolfram.com/language/ref/DistributionFitTest.en.md) * [`EmpiricalDistribution`](https://reference.wolfram.com/language/ref/EmpiricalDistribution.en.md) * [`SmoothKernelDistribution`](https://reference.wolfram.com/language/ref/SmoothKernelDistribution.en.md) * [`QuantityDistribution`](https://reference.wolfram.com/language/ref/QuantityDistribution.en.md) ## Related Guides * [Probability & Statistics with Quantities](https://reference.wolfram.com/language/guide/ProbabilityWithQuantities.en.md) * [`Reliability`](https://reference.wolfram.com/language/guide/Reliability.en.md) * [Probability & Statistics](https://reference.wolfram.com/language/guide/ProbabilityAndStatistics.en.md) * [Random Variables](https://reference.wolfram.com/language/guide/RandomVariables.en.md) * [Survival Analysis](https://reference.wolfram.com/language/guide/SurvivalAnalysis.en.md) * [Statistical Data Analysis](https://reference.wolfram.com/language/guide/Statistics.en.md) * [`Optimization`](https://reference.wolfram.com/language/guide/Optimization.en.md) * [Scientific Data Analysis](https://reference.wolfram.com/language/guide/ScientificDataAnalysis.en.md) * [Robust Descriptive Statistics](https://reference.wolfram.com/language/guide/RobustDescriptiveStatistics.en.md) * [Descriptive Statistics](https://reference.wolfram.com/language/guide/DescriptiveStatistics.en.md) * [Tabular Modeling](https://reference.wolfram.com/language/guide/TabularModeling.en.md) * [Time Series Processing](https://reference.wolfram.com/language/guide/TimeSeries.en.md) * [Actuarial Computation](https://reference.wolfram.com/language/guide/ActuarialComputation.en.md) ## History * [Introduced in 2010 (8.0)](https://reference.wolfram.com/language/guide/SummaryOfNewFeaturesIn80.en.md)