Home Forums Digamma function’s exception, why?

Viewing 1 post (of 1 total)
  • Author
    Posts
  • #1871
    hessinemaaoui
    Keymaster

    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 [2];
        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.