How to model the parameters of AssymetricLaplace and Pareto distribution

Questions version agnostic
1

I have two separate datasets. I want to build a posterior distribution based on Assymetric Laplace for the first and a Pareto distribution for the second.

In both cases, I would take those parameters and then empirically find a distribution that fits each parameter set best. My problem is, how do I go from raw data to getting the parameters needed. Namely, b, kappa and mu for Assymetric Laplace and alpha and m for Pareto?

For reference, if I was modeling a beta distribution, there is an explicit formula to go from mu and sigma of the dataset to alpha and beta. I donā€™t see an analogue for the distributions described above.

As per the request in the comments, I am adding example code for clarity:

from fitter import Fitter
import numpy as np
import pandas as pd
import pymc3 as pm

pymc3_avail_distr = ['normal', 'halfnorm', 'beta', 'expon', 'laplace', 
                    'laplace_asymmetric', 't', 'cauchy', 'halfcauchy',
                    'gamma', 'invgamma', 'lognorm', 'chisquare', 
                    'wald', 'pareto', 'vonmises', 'triangular',
                    'rice', 'logistic']

def get_best_distr(data):
    f = Fitter(data, 
        distributions=pymc3_avail_distr)
    f.fit()
    res = f.get_best()
    return res

df_for_prior = pd.DataFrame(np.random.uniform(low=1700000, high=1900000, size = (22,5)), index = np.arange(1998, 2020))

# transform the values in prior to get the Assymetric LaPlace parameters
# this is the step I am missing

# run fitter function on each parameter to get the best-fitting distribution along with the parameter values to use for specifying the priors

obs = np.array([1.726567e+06, 1.589836e+06, 1.643981e+06, 1.584314e+06])

with pm.Model() as model1:
    b = # some distribution based on results above
    kappa = # 
    mu = #
    
    # Based on distr of raw data
    pred = pm.AssymetricLaplace("modeling_var", b = b, kappa = kappa, mu = mu, observed=obs)

    trace = pm.sample(1000, tune=800)
1 Like
2

Could you elaborate a bit? Do you already have these two models or are you asking about how to build them?

What are ā€œthose parametersā€?

In cases where are no (easy) analytic answers, one tends to builds a model and use an MCMC-based sampling procedure to approximate a posterior which reflects credible values for all your parameters. Thatā€™s part of why I asked my first question above.

That link explains two different specifications of a single distribution. The \mu and \sigma are not ā€œof the datasetā€ any more than the \alpha and \beta are ā€œof the datasetā€.

3

ā€œThose parametersā€ are b, kappa and mu for Assymetric Laplace and alpha and m for Pareto.

Iā€™m asking about how to build these two models. Namely, I specify Assymetric Laplace and Pareto for the priors and want to understand how to specify the priors.

I have a dataset from which I can derive the parameters. If I was specifying my posterior to be a beta distribution, then I cold calculate the alpha and beta parameters from the mean and sigma of the data. Then I can use a library like Fitter to find the distribution that fits that data best and get the relevant parameters for that distribution. However, for the distributions I mention, I donā€™t understand how to calculate the parameters.

That link was intended to show that itā€™s easy to specify Ī± and Ī² in terms of Ī¼ and Ļƒ for the beta distribution.

1 Like
4

Can we see the model so far? Maybe that would help to clarify things.

5

Added an example in the original post. For clarity, I can get a parameter from each column in the prior dataframe. For example mu and sigma. Then I would have 5 mus and sigmas (or whatever other parameters) that I can model distributions from.

6

It sounds like you are trying to do optimization rather than Bayesian inference. If so, then I would suggest using scipyā€™s built-in optimization routines.

Alternatively, it seems like the kind of thing one could do would be the following:

obs = np.array([1.726567e+06, 1.589836e+06, 1.643981e+06, 1.584314e+06])

with pm.Model() as model1:
    kappa = pm.HalfNormal("kappa", 10)
    mu = pm.Normal("mu", mu=np.mean(obs), sigma=np.std(obs))
    b = pm.HalfNormal("b", 3 * np.std(obs))
    
    # Based on distr of raw data
    likelihood = pm.AsymetricLaplace("likelihood", kappa = kappa, mu = mu, b = b, observed=obs)

    idata = pm.sample()
7

Iā€™m not trying to do optimization. Iā€™m trying to get an updated distribution of observed values of which I donā€™t have that many by incorporating historical values (that are a proxy for the measured data) which I do have a lot of.

How do you know that kappa and b follow HalfNormal with the parameters you specified? Also, the priors should be based on the data in df_for_prior, not the observed.

8

If that is the case, then you need to decide how to model your data. If you think your data is distributed as a asymmetric Laplacian, then you need to decide what values of the associated parameters (kappa, mu, b) might be before seeing any data. In the example I presented above, I treated each as uncertain and modeled them as following their own (prior) distribution.

I donā€™t. But kappa and b must take on strictly positive values and using a half normal prior ensures that only positive values are considered. If you believe an alternative is more appropriate, then you should you that.

Below would be a quick way to do what (I think) you are asking for. Itā€™s not clear to me why there are 5 columns or how they differ. So here, Iā€™ll just use 1.

df_for_prior = pd.DataFrame(np.random.uniform(low=1700000, high=1900000, size = (22,5)), index = np.arange(1998, 2020))

obs = np.array([1.726567e+06, 1.589836e+06, 1.643981e+06, 1.584314e+06])

with pm.Model() as model1:
    kappa = pm.HalfNormal("kappa", 10)
    mu = pm.Normal("mu", mu=np.mean(obs), sigma=np.std(obs))
    b = pm.HalfNormal("b", 3 * np.std(obs))
    
    # Based on distr of raw data
    likelihood = pm.AsymetricLaplace("likelihood",
                                     kappa = kappa,
                                     mu = mu,
                                     b = b,
                                     observed=df_for_prior.iloc[:,0])
    likelihood = pm.AsymetricLaplace("likelihood",
                                     kappa = kappa,
                                     mu = mu,
                                     b = b,
                                     observed=obs)

    idata = pm.sample()

Here I just set the priors, updated based on the df_for_prior, then updated further based on obs. If you really have lots of old data, then I would suggest reading through this notebook.

9

Yes, my data is distributed as Assymetric Laplace. The crux of my question was how to specify the associated parameters. Itā€™s not good enough to model them as half-normal, because we actually have data (the prior data) and we should check empirically for the best distribution.

Regarding the data, each column is annual measurements from a physical location that is similar to the target location described by the observed data. Thus, itā€™s not the same, but is a good-enough proxy. So, we have annual data from 5 similar sites. From this, we can build out the distribution of the parameters for the priors.

As I described, if I were modeling a beta distribution for the posterior, then I could calculate the alpha and beta parameters.

Namely, each column would have a mu and a standard deviation. I would then use those to calculate alpha and beta according to the formulas here. Then, I would have 5 alphas and betas (in my real dataset I actually have 20 columns) and could fit a distribution.

The problem is that I donā€™t know how to go from myu and sigma of each column to the parameters in Assymteric Laplace.

10

You certainly can do that, but by no means must you do that. Doing that is equivalent to conducting optimization in which you maximize the likelihood of your data. Suggesting that this procedure produces ā€œthe correctā€ answer is simply false. It is arguably a reasonable answer, but it is one of many reasonable answers.

If you are simply looking for the relationship between the asymmetric Laplacianā€™s parameter values and its moments, you should check the relevant wikipedia page.