statdist.c

speech signal process tools

         10.08,
         12.00,
         13.53,
         16.34,
         19.51,
         21.62,
         24.77,
         27.59,
         31.00,
         33.41,
         40.75
       }
     },
     { "df18",
       { 7.02,
         7.91,
         9.39,
         10.86,
         12.86,
         14.44,
         17.34,
         20.60,
         22.76,
         25.99,
         28.87,
         32.35,
         34.80,
         42.31
       }
     },
     { "df19",
       { 7.63,
         8.57,
         10.12,
         11.65,
         13.72,
         15.35,
         18.34,
         21.69,
         23.90,
         27.20,
         30.14,
         33.69,
         36.19,
         43.82
       }
     },
     { "df20",
       { 8.26,
         9.24,
         10.85,
         12.44,
         14.58,
         16.27,
         19.34,
         22.78,
         25.04,
         28.41,
         31.41,
         35.02,
         37.57,
         45.32
       }
     },
     { "df21",
       { 8.90,
         9.92,
         11.59,
         13.24,
         15.44,
         17.18,
         20.34,
         23.86,
         26.17,
         29.62,
         32.67,
         36.34,
         38.93,
         46.80
       }
     },
     { "df22",
       { 9.54,
         10.60,
         12.34,
         14.04,
         16.31,
         18.10,
         21.24,
         24.94,
         27.30,
         30.81,
         33.92,
         37.66,
         40.29,
         48.27
       }
     },
     { "df23",
       { 10.20,
         11.29,
         13.09,
         14.85,
         17.19,
         19.02,
         22.34,
         26.02,
         28.43,
         32.01,
         35.17,
         38.97,
         41.64,
         49.73
       }
     },
     { "df24",
       { 10.86,
         11.99,
         13.85,
         15.66,
         18.06,
         19.94,
         23.34,
         27.10,
         29.55,
         33.20,
         36.42,
         40.27,
         42.98,
         51.18
       }
     },
     { "df25",
       { 11.52,
         12.70,
         14.61,
         16.47,
         18.94,
         20.87,
         24.34,
         28.17,
         30.68,
         34.38,
         37.65,
         41.57,
         44.31,
         52.62
       }
     },
     { "df26",
       { 12.20,
         13.41,
         15.38,
         17.29,
         19.82,
         21.79,
         25.34,
         29.25,
         31.80,
         35.56,
         38.88,
         42.86,
         45.64,
         54.05
       }
     },
     { "df27",
       { 12.88,
         14.12,
         16.15,
         18.11,
         20.70,
         22.72,
         26.34,
         30.32,
         32.91,
         36.74,
         40.11,
         44.14,
         46.96,
         55.48
       }
     },
     { "df28",
       { 13.56,
         14.85,
         16.93,
         18.94,
         21.59,
         23.65,
         27.34,
         31.39,
         34.03,
         37.92,
         41.34,
         45.42,
         48.28,
         56.89
       }
     },
     { "df29",
       { 14.26,
         15.57,
         17.71,
         19.77,
         22.48,
         24.58,
         28.34,
         32.46,
         35.14,
         39.09,
         42.56,
         46.69,
         49.59,
         58.30
       }
     },
     { "df30",
       { 14.95,
         16.31,
         18.49,
         20.60,
         23.36,
         25.51,
         29.34,
         33.53,
         36.25,
         40.26,
         43.77,
         47.96,
         50.89,
         59.70
       }
     }
    }

};

/**********************************************************/
/*    print to stdout the chi squared distribution table  */
/**********************************************************/
void dump_X2_table(void)
{
    int i,j;
    
    printf("\n\n\t\t\t\tX2 DISTRIBUTION TABLE\n\n");
    printf("\t");
    for (i=MIN_X2_PER;i<MAX_X2_PER+1; i++)
        printf("%s\t",X2.per_str[i]);
    printf("\n\n");
    for (i=MIN_DF; i<MAX_DF+1; i++){
        printf("%s\t",X2.df[i].str);
        for (j=MIN_X2_PER;j<MAX_X2_PER+1; j++)
            printf("%2.3f\t",X2.df[i].level[j]);
        printf("\n");
    }
    printf("\n\n");
}

/**********************************************************************/
/*  given a list of integers, calculate the mean, variance,           */
/*  standard deviation and Z_statistic                                */
/*  Sep 17 1992: changed the variance divisor to be n-1               */
void calc_mean_var_std_dev_Zstat(int *Z_list, int num_Z, float *mean, float *variance, float *std_dev, float *Z_stat)
{
    int i;
    float tmp = 0.0;

    for (i=0;i<num_Z;i++){
       tmp+=Z_list[i];
    }
    *mean = tmp/(float)num_Z;

    tmp = 0.0;
    for (i=0;i<num_Z;i++)
       tmp+= (Z_list[i] - *mean) *  (Z_list[i] - *mean);

    if (tmp != 0.0){
      *variance = tmp / (float)(num_Z - 1);
       *std_dev = sqrt(*variance);
       *Z_stat = *mean / (sqrt(*variance) / sqrt((double)num_Z));
    }
    else
       *variance = *std_dev = *Z_stat = 0.0;
}

/**********************************************************************/
/*  given a list of floats, calculate the mean, variance,             */
/*  standard deviation and Z_statistic                                */
void calc_mean_var_std_dev_Zstat_float(float *Z_list, int num_Z, float *mean, float *variance, float *std_dev, float *Z_stat)
{
    int i;
    float tmp = 0.0;

    for (i=0;i<num_Z;i++){
       tmp+=Z_list[i];
    }
    *mean = tmp/(float)num_Z;

    tmp = 0.0;
    for (i=0;i<num_Z;i++)
       tmp+= (Z_list[i] - *mean) *  (Z_list[i] - *mean);

    if (tmp != 0.0){
       *variance = tmp / (float)(num_Z - 1);
       *std_dev = sqrt(*variance);
       *Z_stat = *mean / (sqrt(*variance) / sqrt((double)num_Z));
    }
    else
       *variance = *std_dev = *Z_stat = 0.0;
}

/**********************************************************************/
/*     a general analysis routine to test whether or not the the      */
/*     Z_statistic indicates a significant difference                 */
int print_Z_analysis(float Z_stat)
{
    int i;

    printf("%s                                  Reject if\n","");
    printf("%s                         Z > (+%1.3f) or Z < (-%1.3f)\n","",
                                                     Z2tail[GEN_Z_PER].z,
                                                     Z2tail[GEN_Z_PER].z);
    printf("\n");
    printf("%s                                  Z = %1.3f\n","",Z_stat);

    printf("\n\n%s\t\tSUMMARY:\n\t\t-------\n\n","");
    if (fabs(Z_stat) > Z2tail[GEN_Z_PER].z){
     printf("\tThere is a significant difference between the test results\n");
        printf("\tusing a %s confidence level.\n",
                                   Z2tail[GEN_Z_PER].str);
        printf("\n");
        printf("\tFurther, the probablity of there being a difference is\n");
        for (i=GEN_Z_PER;i>0;i--)
            if (fabs(Z_stat) < Z2tail[i-1].z)
                break;
        if (i==MAX_Z_PER)
            printf("\tgreater that %s.\n",Z2tail[0].str);
        else
            printf("\tbetween %s to %s.",Z2tail[i].str,Z2tail[i-1].str);
        return(TEST_DIFF);
    }
    else{
        printf("\tThere is no significant difference between the test\n");
        printf("\tresults using a %s confidence level.\n",
                                            Z2tail[GEN_Z_PER].str);
        printf("\n");
        printf("\tFurther, the probablity of there being a difference is\n");
        for (i=GEN_Z_PER;i<MIN_Z_PER;i++)
            if (fabs(Z_stat) > Z2tail[i+1].z)
                break;
        if (i==MIN_Z_PER)
            printf("\tless than %s.\n",Z2tail[i].str);
        else
            printf("\tbetween %s to %s.",Z2tail[i+1].str,Z2tail[i].str);
        return(NO_DIFF);
    }
}

int Z_pass(float Z_stat)
{
    if (fabs(Z_stat) > Z2tail[GEN_Z_PER].z)
        return(TEST_DIFF);
    else
        return(NO_DIFF);
}

/**********************************************************************/
/* calc_two_sample_z_test_float by Brett 5/10/93                      */ 
/*     a general routine to perform a two-sample z test on two lists  */
/*     of floats                                                      */
/**********************************************************************/
void calc_two_sample_z_test_float(float *l1, float *l2, int num_l1, int num_l2, float *Z)
{
  float mean_l1,var_l1,sd_l1,Z_stat_l1;
  float mean_l2,var_l2,sd_l2,Z_stat_l2;
  calc_mean_var_std_dev_Zstat_float(l1,num_l1,&mean_l1,&var_l1,&sd_l1,&Z_stat_l1);
  calc_mean_var_std_dev_Zstat_float(l2,num_l2,&mean_l2,&var_l2,&sd_l2,&Z_stat_l2);
  *Z=((mean_l1-mean_l2)/sqrt(((sd_l1*sd_l1)/(float)num_l1)+((sd_l2*sd_l2)/(float)num_l2)));
  /* print_Z_analysis(*Z);*/
}

/******************************************************************/
/*   Compute the accumulated binomial distribution for the given  */
/*   parameters.  The Formula was taken from                      */
/*   Statistics: Probability, Inference and Decision              */
/*   By Robert L. Winkler and William L. Hays,  Page 206          */
/*                                                                */
/*                        / n \      r      n-r                   */
/*     P( R = r | n, p) = |   |  *  p   *  q                      */
/*                        \ r /                                   */
/*                                                                */
/*        n = number of trials                                    */
/*        R = number of successes in 'n' trials                   */
/*        p = Probability of success                              */
/*        q = 1 - p                                               */
/******************************************************************/
double compute_acc_binomial(int R, int n, double p)
{
    int r, dbg=0;
    double sum=0.0;

    /* printf("Binomial 7,5,0.5 = %f\n",(float)compute_acc_binomial(1,1,0.1));*/

    for (r=0; r <= R; r++) {
        if (dbg) printf("Binomial nCr(%d,%d) = %f / ",
                        n,r,n_CHOOSE_r(n,r));
        if (dbg) printf("[pow(%5.3f,%2d) = %f",
                        p,r,pow(p,(double)r));
        if (dbg) printf(" * pow(%5.3f,%2d) = %f",
                        p,n-r,pow(1.0-p,(double)(n-r)));
        sum += n_CHOOSE_r(n,r) * pow(p,(double)r) * pow(1.0-p,(double)(n-r));	
        if (dbg) printf("] = %f\n",sum);
    }
    return(sum);
}


/*******************************************************************/
/*  Return the result of multipying all the numbers in the sequence*/
/*  from f to to                                                   */
/*******************************************************************/
double seq_mult(int f, int t)
{  
    double sum=1.0;
    int i;
    for (i=f; i<=t; i++)
        sum *= i;
    return(sum);
}

/*******************************************************************/
/* return the result for the N choose R counting Function          */
/*******************************************************************/
double n_CHOOSE_r(int n, int r)
{
    return(seq_mult(r+1,n) / seq_mult(1,n-r));
}