Fixed bug #11227 (T-Test not working with alpha independent variable )
[pspp-builds.git] / src / examine.q
index f8d8ad36ea8821c281e8cf9d6e38a1117309d3d6..d426ecca1e549c5a83dd4d512bfd9b73bcbf0ef9 100644 (file)
@@ -41,6 +41,8 @@ Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA
 #include "hash.h"
 #include "casefile.h"
 #include "factor_stats.h"
+#include "moments.h"
+
 /* (headers) */
 #include "chart.h"
 
@@ -92,48 +94,6 @@ static struct factor *factors=0;
 
 static struct metrics *totals=0;
 
-void
-print_factors(void)
-{
-  struct factor *f = factors;
-
-  while (f) 
-    {
-      struct  factor_statistics **fs = f->fs;
-
-      printf("Factor: %s BY %s\n", 
-            var_to_string(f->indep_var[0]),
-            var_to_string(f->indep_var[1]) );
-
-
-      printf("Contains %d entries\n", hsh_count(f->fstats));
-
-      
-      while (*fs) 
-       {
-         printf("Factor %g; %g\n", (*fs)->id[0].f, (*fs)->id[1].f);
-         
-         /* 
-            printf("Factor %s; %s\n",
-            value_to_string(&(*fs)->id[0], f->indep_var[0]),
-            value_to_string(&(*fs)->id[1], f->indep_var[1]));
-         */
-
-                
-         printf("Sum is %g; ",(*fs)->m[0].sum);
-         printf("N is %g; ",(*fs)->m[0].n);
-         printf("Mean is %g\n",(*fs)->m[0].mean);
-
-         fs++ ;
-       }
-
-      f = f->next;
-    }
-
-  
-}
-
-
 /* Parse the clause specifying the factors */
 static int examine_parse_independent_vars(struct cmd_examine *cmd);
 
@@ -157,7 +117,6 @@ void np_plot(const struct metrics *m, const char *factorname);
 
 
 
-
 /* Per Split function */
 static void run_examine(const struct casefile *cf, void *cmd_);
 
@@ -224,13 +183,29 @@ output_examine(void)
 
       if ( cmd.sbc_plot) 
        {
+         int v;
          if ( cmd.a_plot[XMN_PLT_NPPLOT] ) 
            {
-             int v;
+             for ( v = 0 ; v < n_dependent_vars; ++v ) 
+               np_plot(&totals[v], var_to_string(dependent_vars[v]));
+           }
 
+         if ( cmd.a_plot[XMN_PLT_HISTOGRAM] ) 
+           {
              for ( v = 0 ; v < n_dependent_vars; ++v ) 
-                 np_plot(&totals[v], var_to_string(dependent_vars[v]));
+               {
+                 struct normal_curve normal;
+
+                 normal.N      = totals[v].n;
+                 normal.mean   = totals[v].mean;
+                 normal.stddev = totals[v].stddev;
+                 
+                 histogram_plot(totals[v].histogram, 
+                                var_to_string(dependent_vars[v]),
+                                &normal, 0);
+               }
            }
+
        }
 
 
@@ -254,47 +229,59 @@ output_examine(void)
 
       if ( cmd.sbc_plot) 
        {
-         if ( cmd.a_plot[XMN_PLT_NPPLOT] ) 
+         int v;
+
+         struct factor_statistics **fs = fctr->fs ;
+
+         for ( v = 0 ; v < n_dependent_vars; ++v )
            {
-             int v;
-             for ( v = 0 ; v < n_dependent_vars; ++ v)
+
+             for ( fs = fctr->fs ; *fs ; ++fs ) 
                {
-                 
-                 struct factor_statistics **fs = fctr->fs ;
-                 for ( fs = fctr->fs ; *fs ; ++fs ) 
-                   {
-                     char buf1[100];
-                     char buf2[100];
-                     sprintf(buf1, "%s (",
-                             var_to_string(dependent_vars[v]));
+                 char buf1[100];
+                 char buf2[100];
+                 sprintf(buf1, "%s (",
+                         var_to_string(dependent_vars[v]));
                      
-                     sprintf(buf2, "%s = %s",
-                            var_to_string(fctr->indep_var[0]),
-                            value_to_string(&(*fs)->id[0],fctr->indep_var[0]));
+                 snprintf(buf2, 100, "%s = %s",
+                          var_to_string(fctr->indep_var[0]),
+                          value_to_string(&(*fs)->id[0],fctr->indep_var[0]));
                      
+                 strcat(buf1, buf2);
+                     
+                 if ( fctr->indep_var[1] ) 
+                   {
+                     sprintf(buf2, "; %s = %s)",
+                             var_to_string(fctr->indep_var[1]),
+                             value_to_string(&(*fs)->id[1],
+                                             fctr->indep_var[1]));
                      strcat(buf1, buf2);
+                   }
+                 else
+                   {
+                     strcat(buf1, ")");
+                   }
 
-                     
-                     if ( fctr->indep_var[1] ) 
-                       {
-                         sprintf(buf2, "; %s = %s)",
-                                 var_to_string(fctr->indep_var[1]),
-                                 value_to_string(&(*fs)->id[1],
-                                                 fctr->indep_var[1]));
-                         strcat(buf1, buf2);
-                       }
-                     else
-                       {
-                         strcat(buf1, ")");
-                       }
-
-                     np_plot(&(*fs)->m[v],buf1);
+                 if ( cmd.a_plot[XMN_PLT_NPPLOT] ) 
+                   np_plot(&(*fs)->m[v],buf1);
 
+                 
+                 if ( cmd.a_plot[XMN_PLT_HISTOGRAM] ) 
+                   {
+                     struct normal_curve normal;
+
+                     normal.N      = (*fs)->m[v].n;
+                     normal.mean   = (*fs)->m[v].mean;
+                     normal.stddev = (*fs)->m[v].stddev;
+                 
+                     histogram_plot((*fs)->m[v].histogram, 
+                                    buf1,  &normal, 0);
                    }
                  
-               }
+               } /* for ( fs .... */
+
+           } /* for ( v = 0 ..... */
 
-           }
        }
 
       fctr = fctr->next;
@@ -1289,7 +1276,7 @@ populate_descriptives(struct tab_table *tbl, int col, int row,
   tab_float (tbl, col + 3,
             row,
             TAB_CENTER,
-            m->stderr,
+            m->se_mean,
             8,3);
   
 
@@ -1307,7 +1294,7 @@ populate_descriptives(struct tab_table *tbl, int col, int row,
   tab_float (tbl, col + 2,
             row + 1,
             TAB_CENTER,
-            m->mean - t * m->stderr
+            m->mean - t * m->se_mean
             8,3);
 
   tab_text (tbl, col + 1,  
@@ -1319,7 +1306,7 @@ populate_descriptives(struct tab_table *tbl, int col, int row,
   tab_float (tbl, col + 2,
             row + 2,
             TAB_CENTER,
-            m->mean + t * m->stderr
+            m->mean + t * m->se_mean
             8,3);
 
   tab_text (tbl, col, 
@@ -1408,14 +1395,42 @@ populate_descriptives(struct tab_table *tbl, int col, int row,
            TAB_LEFT | TAT_TITLE,
            _("Skewness"));
 
+
+  tab_float (tbl, col + 2,
+            row + 11,
+            TAB_CENTER,
+            m->skewness,
+            8,3);
+
+  /* stderr of skewness */
+  tab_float (tbl, col + 3,
+            row + 11,
+            TAB_CENTER,
+            calc_seskew(m->n),
+            8,3);
+
+
   tab_text (tbl, col, 
            row + 12,
            TAB_LEFT | TAT_TITLE,
            _("Kurtosis"));
-}
 
 
+  tab_float (tbl, col + 2,
+            row + 12,
+            TAB_CENTER,
+            m->kurtosis,
+            8,3);
+
+  /* stderr of kurtosis */
+  tab_float (tbl, col + 3,
+            row + 12,
+            TAB_CENTER,
+            calc_sekurt(m->n),
+            8,3);
+
 
+}
 
 
 /* Plot the normal and detrended normal plots for m
@@ -1432,9 +1447,6 @@ np_plot(const struct metrics *m, const char *factorname)
   /* Detrended Normal Plot */
   struct chart dnp_chart;
 
-  const struct weighted_value *wv = *(m->wvp);
-
-
   /* The slope and intercept of the ideal normal probability line */
   const double slope = 1.0 / m->stddev;
   const double intercept = - m->mean / m->stddev;
@@ -1454,8 +1466,9 @@ np_plot(const struct metrics *m, const char *factorname)
   chart_write_xlabel(&dnp_chart, _("Observed Value"));
   chart_write_ylabel(&dnp_chart, _("Dev from Normal"));
 
-  yfirst = gsl_cdf_ugaussian_Pinv (wv[0].rank / ( m->n + 1));
-  ylast =  gsl_cdf_ugaussian_Pinv (wv[m->n_data-1].rank / ( m->n + 1));
+  yfirst = gsl_cdf_ugaussian_Pinv (m->wvp[0]->rank / ( m->n + 1));
+  ylast =  gsl_cdf_ugaussian_Pinv (m->wvp[m->n_data-1]->rank / ( m->n + 1));
+
 
   {
     /* Need to make sure that both the scatter plot and the ideal fit into the
@@ -1464,17 +1477,13 @@ np_plot(const struct metrics *m, const char *factorname)
     double x_upper = max(m->max, (ylast  - intercept) / slope) ;
     double slack = (x_upper - x_lower)  * 0.05 ;
 
-    chart_write_xscale(&np_chart, x_lower  - slack, x_upper + slack,
-                      chart_rounded_tick((m->max - m->min) / 5.0));
+    chart_write_xscale(&np_chart, x_lower - slack, x_upper + slack, 5);
 
-
-    chart_write_xscale(&dnp_chart, m->min, m->max,
-                      chart_rounded_tick((m->max - m->min) / 5.0));
+    chart_write_xscale(&dnp_chart, m->min, m->max, 5);
 
   }
 
-  chart_write_yscale(&np_chart, yfirst, ylast, 
-                    chart_rounded_tick((ylast - yfirst)/5.0) );
+  chart_write_yscale(&np_chart, yfirst, ylast, 5);
 
   {
   /* We have to cache the detrended data, beacause we need to 
@@ -1485,21 +1494,19 @@ np_plot(const struct metrics *m, const char *factorname)
   double d_min = DBL_MAX;
   for ( i = 0 ; i < m->n_data; ++i ) 
     {
-      const double ns = gsl_cdf_ugaussian_Pinv (wv[i].rank / ( m->n + 1));
+      const double ns = gsl_cdf_ugaussian_Pinv (m->wvp[i]->rank / ( m->n + 1));
 
-      chart_datum(&np_chart, 0, wv[i].v.f, ns);
+      chart_datum(&np_chart, 0, m->wvp[i]->v.f, ns);
 
-      d_data[i] = (wv[i].v.f - m->mean) / m->stddev  - ns;
+      d_data[i] = (m->wvp[i]->v.f - m->mean) / m->stddev  - ns;
    
       if ( d_data[i] < d_min ) d_min = d_data[i];
       if ( d_data[i] > d_max ) d_max = d_data[i];
     }
-
-  chart_write_yscale(&dnp_chart, d_min, d_max, 
-                    chart_rounded_tick((d_max - d_min) / 5.0));
+  chart_write_yscale(&dnp_chart, d_min, d_max, 5);
 
   for ( i = 0 ; i < m->n_data; ++i ) 
-      chart_datum(&dnp_chart, 0, wv[i].v.f, d_data[i]);
+      chart_datum(&dnp_chart, 0, m->wvp[i]->v.f, d_data[i]);
 
   free(d_data);
   }