Home Forums Digamma function’s exception, why?

• This topic has 1 reply, 2 voices, and was last updated 10 years, 5 months ago by hessinemaaoui.
Viewing 1 post (of 1 total)
• Author
Posts
• #1871

Here is the code:

import com.numericalmethod.suanshu.analysis.uniroot.Newton;
import com.numericalmethod.suanshu.analysis.function.rn2r1.UnivariateRealFunction;
import com.numericalmethod.suanshu.analysis.function.special.Digamma;

public class Trial6{

public static void main(String[] args) {
final double [] v = new double []{0.52392859 ,0.74728222, 0.35648189, 1.50509941, 0.09642664, 2.60895326,  0.12903240, 1.58223482, 0.08668512, 0.47263063, 0.51582472, 0.03652849, 2.33088775, 2.56285203, 0.66409021, 0.10836313, 2.11359703, 0.73836191, 0.28691570, 0.19404304, 0.22538941, 0.04399578, 1.26417153, 3.86406407, 0.12689116, 0.16198531, 0.94862916, 1.10415084, 3.39500931, 0.62369945, 3.70865055, 1.13927840, 0.12278842, 1.17936096, 2.90563559, 1.23769602};
double [] a = new double ;
a=gamma_mle(v);
Digamma digamma = new Digamma();
}

public static double [] gamma_mle (final double []y) {

UnivariateRealFunction fn = new UnivariateRealFunction() {

@Override
public double evaluate(double x) {
Digamma digamma = new Digamma();
return Math.log(x)-digamma.evaluate(x)-Math.log(sum(y)/y.length+sum(log(y))/y.length);
}
};
Newton newton = new Newton(fn, 1e-5);
double alpha = newton.solve(100, 1);
double lambda= alpha*y.length/sum(y);
double [] vectorofpara = new double []{alpha,lambda};
return vectorofpara;

}

public static double sum (double[] x) {
int i,j;
double sum=0;

for ( j = 0; j < x.length; j++) {
sum=sum+x[j];
}
return sum ;

public static double []  log(double[] x) {
int j;
for ( j = 0; j < x.length; j++) {
x[j]= Math.log(x[j]);
}
return x;

}
}

The program is implemented to test the maximum likelihood estimate of gamma distribution (alpha,lambda) (supposed to be 1) , given some observations (v) generated by other methods from the distribution gamma(1,1), the algorithm for finding the mle of the parameters is based on the following link from wiki : http://en.wikipedia.org/wiki/Gamma_distribution , but an error occurs due to exceptions in digamma.evaluate, how should I amend the program?

Viewing 1 post (of 1 total)
• You must be logged in to reply to this topic.