add NIST anova tests

PDLPorters · Jan 27, 2025 · f4d6499 · f4d6499
1 parent fc46091
commit f4d6499
Showing 1 changed file with 67 additions and 0 deletions.
diff --git a/t/glm.t b/t/glm.t
@@ -295,6 +295,73 @@ is_pdl $b_bad->dvrs(ones(6) * .5), pdl( 'BAD -1.17741002251547 -1.17741002251547
   });
 }
 
+{ # anova_nist_low
+  # data from https://www.itl.nist.gov/div898/strd/anova/SmLs01.html
+  #  1   2   3   4   5   6   7   8   9
+  # 1.4 1.3 1.5 1.3 1.5 1.3 1.5 1.3 1.5
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # 1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4
+  # 1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6
+  # Certified Values:
+  # Source of                  Sums of               Mean
+  # Variation          df      Squares              Squares             F Statistic
+  # Between Treatment   8 1.68000000000000E+00 2.10000000000000E-01 2.10000000000000E+01
+  # Within Treatment  180 1.80000000000000E+00 1.00000000000000E-02
+  #  Certified R-Squared 4.82758620689655E-01
+  #  Certified Residual
+  #  Standard Deviation  1.00000000000000E-01
+  my $data = pdl '[
+   [1.4 1.3 1.5 1.3 1.5 1.3 1.5 1.3 1.5]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+   [1.3 1.2 1.4 1.2 1.4 1.2 1.4 1.2 1.4]
+   [1.5 1.4 1.6 1.4 1.6 1.4 1.6 1.4 1.6]
+  ]';
+  my %m = $data->flat->anova(sequence(9)->dummy(1,21)->flat);
+  test_stats_cmp(\%m, {
+    '| IV_0 | ms' => 0.21,
+    '| IV_0 | ss' => 1.68,
+    F => 21,
+    F_df => pdl(8, 180),
+    ss_model => 1.68,
+    ss_residual => 1.8,
+    ms_residual => 0.01,
+  });
+}
+
 is_pdl pdl('BAD 1 2 3 4; BAD BAD BAD BAD BAD')->fill_m, pdl('2.5 1 2 3 4; 0 0 0 0 0'), 'fill_m nan to bad';
 is_pdl pdl([1,1,1], [2,2,2])->stddz, zeroes(3,2), 'stddz nan vs bad';