Bayesian update: Weibull dist. parameters

UQ with UQLab Community Q&A and How To
1

Hi all,

Assume X\sim Weibull(\alpha,\beta) with observations \mathcal{X}=\left\{x_{1},x_{2},\cdots\right\}. Prior distributions for \alpha and \beta are as follows:

\pi(\alpha)=\mathcal{N}(13.6, 1.36)
\pi(\beta)=\mathcal{N}(10.2, 1.02)

My goal is to update the values of \alpha and \beta based on the observations \mathcal{X}.

Dear @paulremo, you discussed the solution and limitations for a simple case in the topic Bayesian inversion using hierarchical prior distributions in the UQlab. What simplifications and/or modifications make it possible to update Weibull distribution parameters in the current version of UQLab?

Any feedback would be greatly appreciated.

Thanks for your consideration.

2

Hi @Soheil_Radfar

The Bayesian inversion module by default assumes that the data are realizations of random variables distributed according to normal distributions. In your case, where \mathcal{X} are realizations of X\sim\text{Weibull}(\alpha,\beta), you need to specify a user-defined likelihood function (have a look at the Bayesian inversion user manual or Example 6 of the Bayesian inversion module).

You already know the priors for \alpha and \beta and the likelihood function is just a product of Weibull PDFs at the data \mathcal{X}.

Let me know if you need further help!

3

Hi @paulremo
Thank you very much for your helpful explanation.

According to your comment, I wrote the custom likelihood function as follows:

function logL = uq_WeibullLogLik(params, y)

measurements = y;

nReal = size(params,1); % number of queried realizations

Ny = size(measurements);
logL = zeros(nReal,1);

    for ii = 1:nReal
    logLikeli = 0;
    
        for jj = 1:Ny
            logLikeli = logLikeli + log(wblpdf(measurements(jj,:),...
                modelParams(ii,1),modelParams(ii,2)));
        end
    
    logL(ii) = logLikeli ;
    end
end

Is it implemented correctly?

4

Hi @Soheil_Radfar

Looks good :blush:! But I think modelparams should be params and for efficiency you should avoid loops in Matlab whenever possible and write the function in vector form instead. I think the following should work:

function logL = uq_WeibullLogLik(params, y)

measurements = y;

nReal = size(params,1); % number of queried realizations

Ny = size(measurements, 2); % should be 1xNy, please double check
logL = zeros(nReal,1);

for ii = 1:nReal
    logLikeli = log(wblpdf(measurements, modelParams(ii,1), modelParams(ii,2)));
        
    logL(ii) = sum(logLikeli);
end
5

Awesome!
Thank you very much @paulremo, your guidance is greatly appreciated.

Here is the final version with some minor modifications:

function logL = uq_WeibullLogLik(params, y)

measurements = y;

nReal = size(params,1); % number of queried realizations

logL = zeros(nReal,1);

for ii = 1:nReal
    logLikeli = log(wblpdf(measurements, params(ii,1), params(ii,2)));
        
    logL(ii) = sum(logLikeli);
end
1 Like