diff --git a/dataframe-core/src-internal/DataFrame/Internal/AggKernel.hs b/dataframe-core/src-internal/DataFrame/Internal/AggKernel.hs
index 185df392..a01ef9f0 100644
--- a/dataframe-core/src-internal/DataFrame/Internal/AggKernel.hs
+++ b/dataframe-core/src-internal/DataFrame/Internal/AggKernel.hs
@@ -217,7 +217,8 @@ varScatter takeSqrt g nGroups v = runST $ do
             | otherwise = do
                 c <- VUM.unsafeRead cnt k
                 mm <- VUM.unsafeRead m2 k
-                let var = if c < 2 then 0 else mm / fromIntegral (c - 1)
+                -- NaN at n = 1, matching computeVariance
+                let var = if c < 2 then 0 / 0 else mm / fromIntegral (c - 1)
                 VUM.unsafeWrite out k (if takeSqrt then sqrt var else var)
                 fin (k + 1)
     fin 0
diff --git a/dataframe-core/src-internal/DataFrame/Internal/AggKernelPar.hs b/dataframe-core/src-internal/DataFrame/Internal/AggKernelPar.hs
index db8ae7a7..cd10c924 100644
--- a/dataframe-core/src-internal/DataFrame/Internal/AggKernelPar.hs
+++ b/dataframe-core/src-internal/DataFrame/Internal/AggKernelPar.hs
@@ -290,7 +290,8 @@ varPar takeSqrt vis offs nGroups v caps bounds = do
             | otherwise = do
                 c <- VUM.unsafeRead cnt k
                 mm <- VUM.unsafeRead m2 k
-                let var = if c < 2 then 0 else mm / fromIntegral (c - 1)
+                -- NaN at n = 1, matching computeVariance
+                let var = if c < 2 then 0 / 0 else mm / fromIntegral (c - 1)
                 VUM.unsafeWrite out k (if takeSqrt then sqrt var else var)
                 fin (k + 1)
     fin 0
diff --git a/dataframe-core/src-internal/DataFrame/Internal/Column.hs b/dataframe-core/src-internal/DataFrame/Internal/Column.hs
index 1395b9b4..19784282 100644
--- a/dataframe-core/src-internal/DataFrame/Internal/Column.hs
+++ b/dataframe-core/src-internal/DataFrame/Internal/Column.hs
@@ -209,6 +209,28 @@ columnBitmap (UnboxedColumn bm _) = bm
 columnBitmap (PackedText bm _) = bm
 columnBitmap (MergedColumn _ _) = Nothing
 
+{- | Drops the null values in a nullable column: those slots hold a sentinel,
+not a value. Identity on columns without a bitmap. 'VG.ifilter' inspects only
+the index, so boxed error thunks at null slots are never forced.
+-}
+dropNulls :: Column -> Column
+dropNulls (BoxedColumn (Just bm) xs) =
+    BoxedColumn Nothing (VG.ifilter (\i _ -> bitmapTestBit bm i) xs)
+dropNulls (UnboxedColumn (Just bm) xs) =
+    UnboxedColumn Nothing (VG.ifilter (\i _ -> bitmapTestBit bm i) xs)
+dropNulls c@(PackedText (Just _) _) = dropNulls (materializePacked c)
+dropNulls c = c
+{-# INLINE dropNulls #-}
+
+{- | 'dropNulls', unless the view type @a@ is @Maybe@-headed: a @Maybe@-typed
+view encodes the nulls as values, so the column passes through untouched.
+-}
+dropNullsExceptMaybe :: forall a. (Typeable a) => Column -> Column
+dropNullsExceptMaybe c = case typeRep @a of
+    App m _ | Just HRefl <- eqTypeRep m (typeRep @Maybe) -> c
+    _ -> dropNulls c
+{-# INLINE dropNullsExceptMaybe #-}
+
 {- | Decode a 'PackedText' into a @BoxedColumn Text@ (bit-identical to
 materializing at freeze). Identity on every other column.
 -}
diff --git a/dataframe-learn/src/DataFrame/Metrics.hs b/dataframe-learn/src/DataFrame/Metrics.hs
index 3982efba..c23a6d21 100644
--- a/dataframe-learn/src/DataFrame/Metrics.hs
+++ b/dataframe-learn/src/DataFrame/Metrics.hs
@@ -1,3 +1,4 @@
+{-# LANGUAGE OverloadedStrings #-}
 {-# LANGUAGE ScopedTypeVariables #-}
 {-# LANGUAGE TypeApplications #-}
 
@@ -45,6 +46,7 @@ import Data.Ord (comparing)
 import qualified Data.Text as T
 import qualified Data.Vector.Unboxed as VU
 
+import DataFrame.Errors (DataFrameException (..))
 import DataFrame.Internal.Column (TypedColumn (..), toVector)
 import DataFrame.Internal.DataFrame (DataFrame)
 import DataFrame.Internal.Expression (Expr)
@@ -73,15 +75,20 @@ columnOf df e = case interpret @Double df e of
     Right (TColumn c) -> fromRight VU.empty (toVector @Double @VU.Vector c)
     Left err -> throw err
 
-n2 :: VU.Vector Double -> Double
-n2 = fromIntegral . VU.length
+{- | Compared pairs: 'VU.zipWith' truncates to the shorter vector, so every
+mean below divides by this, never by the length of 'truth' alone.
+-}
+nCompared :: VU.Vector Double -> VU.Vector Double -> Double
+nCompared preds truth = fromIntegral (min (VU.length preds) (VU.length truth))
 
 -- | Mean squared error.
 mse :: Metric
 mse preds truth
-    | VU.null truth = 0
+    | n == 0 = throw (EmptyDataSetException "mse")
     | otherwise =
-        VU.sum (VU.zipWith (\p t -> (p - t) ^ (2 :: Int)) preds truth) / n2 truth
+        VU.sum (VU.zipWith (\p t -> (p - t) ^ (2 :: Int)) preds truth) / n
+  where
+    n = nCompared preds truth
 
 -- | Root mean squared error.
 rmse :: Metric
@@ -90,30 +97,37 @@ rmse preds truth = sqrt (mse preds truth)
 -- | Mean absolute error.
 mae :: Metric
 mae preds truth
-    | VU.null truth = 0
-    | otherwise = VU.sum (VU.zipWith (\p t -> abs (p - t)) preds truth) / n2 truth
+    | n == 0 = throw (EmptyDataSetException "mae")
+    | otherwise = VU.sum (VU.zipWith (\p t -> abs (p - t)) preds truth) / n
+  where
+    n = nCompared preds truth
 
 -- | Coefficient of determination @R²@.
 r2 :: Metric
 r2 preds truth
-    | VU.null truth || ssTot == 0 = 0
+    | n == 0 = throw (EmptyDataSetException "r2")
+    | ssTot == 0 = 0
     | otherwise = 1 - ssRes / ssTot
   where
-    mean = VU.sum truth / n2 truth
-    ssRes = VU.sum (VU.zipWith (\p t -> (t - p) ^ (2 :: Int)) preds truth)
-    ssTot = VU.sum (VU.map (\t -> (t - mean) ^ (2 :: Int)) truth)
+    n = nCompared preds truth
+    truth' = VU.take (min (VU.length preds) (VU.length truth)) truth
+    mean = VU.sum truth' / n
+    ssRes = VU.sum (VU.zipWith (\p t -> (t - p) ^ (2 :: Int)) preds truth')
+    ssTot = VU.sum (VU.map (\t -> (t - mean) ^ (2 :: Int)) truth')
 
 -- | Fraction of exact matches.
 accuracy :: Metric
 accuracy preds truth
-    | VU.null truth = 0
+    | n == 0 = throw (EmptyDataSetException "accuracy")
     | otherwise =
-        fromIntegral (VU.length (VU.filter id (VU.zipWith (==) preds truth))) / n2 truth
+        fromIntegral (VU.length (VU.filter id (VU.zipWith (==) preds truth))) / n
+  where
+    n = nCompared preds truth
 
 -- | Binary log loss; probabilities clamped away from @0@/@1@.
 logLoss :: Metric
 logLoss probs truth
-    | VU.null truth = 0
+    | n == 0 = throw (EmptyDataSetException "logLoss")
     | otherwise =
         negate
             ( VU.sum
@@ -123,8 +137,9 @@ logLoss probs truth
                     truth
                 )
             )
-            / n2 truth
+            / n
   where
+    n = nCompared probs truth
     clampP p = max 1e-15 (min (1 - 1e-15) p)
 
 -- | Averaging strategy for multiclass precision/recall/F1.
diff --git a/dataframe-learn/tests-internal/Learn/EdgeCases.hs b/dataframe-learn/tests-internal/Learn/EdgeCases.hs
index 6a991570..a3012c0e 100644
--- a/dataframe-learn/tests-internal/Learn/EdgeCases.hs
+++ b/dataframe-learn/tests-internal/Learn/EdgeCases.hs
@@ -28,7 +28,7 @@ import DataFrame.LinearModel
 import DataFrame.LinearSolver (sigmoid)
 import DataFrame.PCA
 
-import DataFrame.Internal.Statistics (correlation', variance')
+import DataFrame.Internal.Statistics (correlation', meanSquaredError, variance')
 
 import Test.HUnit
 
@@ -124,14 +124,13 @@ testVarianceConstant = TestCase $ do
     let v = variance' (VU.replicate 100 (7.0 :: Double))
     assertEqual "variance of constant column is 0" 0 v
 
-{- Variance of fewer than two samples is defined to be 0 (computeVariance guard),
-   not NaN from a /0. -}
+{- Sample variance of one observation is undefined: NaN, so a singleton
+   group can never look as tight as a genuinely constant column. -}
 testVarianceSingleton :: Test
 testVarianceSingleton = TestCase $ do
-    assertEqual
-        "variance of one sample is 0"
-        0
-        (variance' (VU.fromList [3.5 :: Double]))
+    assertBool
+        "variance of one sample is NaN"
+        (isNaN (variance' (VU.fromList [3.5 :: Double])))
 
 {- Correlation of a perfectly linear pair is exactly +1 (and -1 reversed),
    computed stably. y = 2x+1 over a spread of x. -}
@@ -169,6 +168,27 @@ testCorrelationTooFew = TestCase $ do
         Nothing
         (correlation' (VU.fromList [1]) (VU.fromList [2]))
 
+{- meanSquaredError refuses length mismatches and empty inputs rather than
+   averaging over terms it never summed (or indexing out of bounds). -}
+testMeanSquaredErrorGuards :: Test
+testMeanSquaredErrorGuards = TestCase $ do
+    assertEqual
+        "mse of mismatched lengths is Nothing"
+        Nothing
+        (meanSquaredError (VU.fromList [0, 0, 0, 0]) (VU.fromList [2, 2]))
+    assertEqual
+        "mse with the longer prediction does not index out of bounds"
+        Nothing
+        (meanSquaredError (VU.fromList [1]) (VU.fromList [1, 2, 3]))
+    assertEqual
+        "mse of empty inputs is Nothing"
+        Nothing
+        (meanSquaredError VU.empty VU.empty)
+    assertEqual
+        "mse of equal-length inputs is the plain mean"
+        (Just 4.0)
+        (meanSquaredError (VU.fromList [0, 0]) (VU.fromList [2, 2]))
+
 -- ===========================================================================
 -- Category 8: stability inside the model expr layer
 -- ===========================================================================
@@ -425,6 +445,7 @@ tests =
     , testCorrelationPerfect
     , testCorrelationConstantColumnIsNaN
     , testCorrelationTooFew
+    , testMeanSquaredErrorGuards
     , testLogisticProbsExtremeFeatures
     , testOLSOneRow
     , testLogisticSingleClass
diff --git a/dataframe-operations/src-internal/DataFrame/Internal/Statistics.hs b/dataframe-operations/src-internal/DataFrame/Internal/Statistics.hs
index 0c80ad4a..ef929d83 100644
--- a/dataframe-operations/src-internal/DataFrame/Internal/Statistics.hs
+++ b/dataframe-operations/src-internal/DataFrame/Internal/Statistics.hs
@@ -73,7 +73,9 @@ varianceStep (VarAcc !n !meanVal !m2) !x =
 
 computeVariance :: VarAcc -> Double
 computeVariance (VarAcc !n _ !m2)
-    | n < 2 = 0 -- or error "variance of <2 samples"
+    | n == 0 = throw $ EmptyDataSetException "variance"
+    -- undefined at n = 1: NaN, not 0
+    | n < 2 = 0 / 0
     | otherwise = m2 / fromIntegral (n - 1)
 {-# INLINE computeVariance #-}
 
@@ -106,7 +108,7 @@ skewnessStep (SkewAcc !n !meanVal !m2 !m3) !x' =
 computeSkewness :: SkewAcc -> Double
 computeSkewness (SkewAcc n _ m2 m3)
     | n < 3 = 0 -- or error "skewness of <3 samples"
-    | otherwise = (sqrt (fromIntegral n - 1) * m3) / sqrt (m2 ^ (3 :: Int))
+    | otherwise = (sqrt (fromIntegral n) * m3) / sqrt (m2 ^ (3 :: Int))
 {-# INLINE computeSkewness #-}
 
 skewness' :: (VU.Unbox a, Real a, Num a) => VU.Vector a -> Double
@@ -202,11 +204,14 @@ interQuartileRange' samp =
 {-# INLINE interQuartileRange' #-}
 
 meanSquaredError :: VU.Vector Double -> VU.Vector Double -> Maybe Double
-meanSquaredError target prediction =
-    let
-        squareDiff = VU.ifoldl' (\sq i e -> (e - target VU.! i) ^ (2 :: Int) + sq) 0 prediction
-     in
-        Just $ squareDiff / fromIntegral (max (VU.length target) (VU.length prediction))
+meanSquaredError target prediction
+    | VU.length target /= VU.length prediction = Nothing
+    | VU.null target = Nothing
+    | otherwise =
+        Just
+            ( VU.sum (VU.zipWith (\t p -> (p - t) ^ (2 :: Int)) target prediction)
+                / fromIntegral (VU.length target)
+            )
 {-# INLINE meanSquaredError #-}
 
 mutualInformationBinned ::
diff --git a/dataframe-operations/src/DataFrame/Operations/Statistics.hs b/dataframe-operations/src/DataFrame/Operations/Statistics.hs
index 7322c0dc..1526d500 100644
--- a/dataframe-operations/src/DataFrame/Operations/Statistics.hs
+++ b/dataframe-operations/src/DataFrame/Operations/Statistics.hs
@@ -118,7 +118,7 @@ mean (Col name) df = case _getColumnAsDouble name df of
     Nothing -> error "[INTERNAL ERROR] Column is non-numeric"
 mean expr df = case interpret df expr of
     Left e -> throw e
-    Right (TColumn col) -> case toUnboxedVector @a col of
+    Right (TColumn col) -> case toUnboxedVector @a (dropNulls col) of
         Left e -> throw e
         Right xs -> mean' xs
 
@@ -136,12 +136,15 @@ meanMaybe expr df = case interpret @(Maybe a) df expr of
 -- | Calculates the median of a given column as a standalone value.
 median ::
     forall a. (Columnable a, Real a, VU.Unbox a) => Expr a -> DataFrame -> Double
-median (Col name) df = case columnAsUnboxedVector (Col @a name) df of
-    Right xs -> median' xs
-    Left e -> throw e
+median (Col name) df = case getColumn name df of
+    Just col -> case toUnboxedVector @a (dropNulls col) of
+        Right xs -> median' xs
+        Left e -> throw e
+    Nothing ->
+        throw $ ColumnsNotFoundException [name] "median" (M.keys $ columnIndices df)
 median expr df = case interpret df expr of
     Left e -> throw e
-    Right (TColumn col) -> case toUnboxedVector @a col of
+    Right (TColumn col) -> case toUnboxedVector @a (dropNulls col) of
         Left e -> throw e
         Right xs -> median' xs
 
@@ -161,12 +164,15 @@ medianMaybe expr df = case interpret @(Maybe a) df expr of
 percentile ::
     forall a.
     (Columnable a, Real a, VU.Unbox a) => Int -> Expr a -> DataFrame -> Double
-percentile n (Col name) df = case columnAsUnboxedVector (Col @a name) df of
-    Right xs -> percentile' n xs
-    Left e -> throw e
+percentile n (Col name) df = case getColumn name df of
+    Just col -> case toUnboxedVector @a (dropNulls col) of
+        Right xs -> percentile' n xs
+        Left e -> throw e
+    Nothing ->
+        throw $ ColumnsNotFoundException [name] "percentile" (M.keys $ columnIndices df)
 percentile n expr df = case interpret df expr of
     Left e -> throw e
-    Right (TColumn col) -> case toUnboxedVector @a col of
+    Right (TColumn col) -> case toUnboxedVector @a (dropNulls col) of
         Left e -> throw e
         Right xs -> percentile' n xs
 
@@ -174,36 +180,47 @@ percentile n expr df = case interpret df expr of
 genericPercentile ::
     forall a.
     (Columnable a, Ord a) => Int -> Expr a -> DataFrame -> a
-genericPercentile n (Col name) df = case columnAsVector (Col @a name) df of
-    Right xs -> percentileOrd' n xs
-    Left e -> throw e
+genericPercentile n (Col name) df = case getColumn name df of
+    Just col -> case toVector @a (dropNullsExceptMaybe @a col) of
+        Right xs -> percentileOrd' n xs
+        Left e -> throw e
+    Nothing ->
+        throw $
+            ColumnsNotFoundException [name] "genericPercentile" (M.keys $ columnIndices df)
 genericPercentile n expr df = case interpret df expr of
     Left e -> throw e
-    Right (TColumn col) -> case toVector @a col of
+    Right (TColumn col) -> case toVector @a (dropNullsExceptMaybe @a col) of
         Left e -> throw e
         Right xs -> percentileOrd' n xs
 
 -- | Calculates the standard deviation of a given column as a standalone value.
 standardDeviation ::
     forall a. (Columnable a, Real a, VU.Unbox a) => Expr a -> DataFrame -> Double
-standardDeviation (Col name) df = case columnAsUnboxedVector (Col @a name) df of
-    Right xs -> (sqrt . variance') xs
-    Left e -> throw e
+standardDeviation (Col name) df = case getColumn name df of
+    Just col -> case toUnboxedVector @a (dropNulls col) of
+        Right xs -> (sqrt . variance') xs
+        Left e -> throw e
+    Nothing ->
+        throw $
+            ColumnsNotFoundException [name] "standardDeviation" (M.keys $ columnIndices df)
 standardDeviation expr df = case interpret df expr of
     Left e -> throw e
-    Right (TColumn col) -> case toUnboxedVector @a col of
+    Right (TColumn col) -> case toUnboxedVector @a (dropNulls col) of
         Left e -> throw e
         Right xs -> (sqrt . variance') xs
 
 -- | Calculates the skewness of a given column as a standalone value.
 skewness ::
     forall a. (Columnable a, Real a, VU.Unbox a) => Expr a -> DataFrame -> Double
-skewness (Col name) df = case columnAsUnboxedVector (Col @a name) df of
-    Right xs -> skewness' xs
-    Left e -> throw e
+skewness (Col name) df = case getColumn name df of
+    Just col -> case toUnboxedVector @a (dropNulls col) of
+        Right xs -> skewness' xs
+        Left e -> throw e
+    Nothing ->
+        throw $ ColumnsNotFoundException [name] "skewness" (M.keys $ columnIndices df)
 skewness expr df = case interpret df expr of
     Left e -> throw e
-    Right (TColumn col) -> case toUnboxedVector @a col of
+    Right (TColumn col) -> case toUnboxedVector @a (dropNulls col) of
         Left e -> throw e
         Right xs -> skewness' xs
 
@@ -215,42 +232,51 @@ variance (Col name) df = case _getColumnAsDouble name df of
     Nothing -> error "[INTERNAL ERROR] Column is non-numeric"
 variance expr df = case interpret df expr of
     Left e -> throw e
-    Right (TColumn col) -> case toUnboxedVector @a col of
+    Right (TColumn col) -> case toUnboxedVector @a (dropNulls col) of
         Left e -> throw e
         Right xs -> variance' xs
 
 -- | Calculates the inter-quartile range of a given column as a standalone value.
 interQuartileRange ::
     forall a. (Columnable a, Real a, VU.Unbox a) => Expr a -> DataFrame -> Double
-interQuartileRange (Col name) df = case columnAsUnboxedVector (Col @a name) df of
-    Right xs -> interQuartileRange' xs
-    Left e -> throw e
+interQuartileRange (Col name) df = case getColumn name df of
+    Just col -> case toUnboxedVector @a (dropNulls col) of
+        Right xs -> interQuartileRange' xs
+        Left e -> throw e
+    Nothing ->
+        throw $
+            ColumnsNotFoundException [name] "interQuartileRange" (M.keys $ columnIndices df)
 interQuartileRange expr df = case interpret df expr of
     Left e -> throw e
-    Right (TColumn col) -> case toUnboxedVector @a col of
+    Right (TColumn col) -> case toUnboxedVector @a (dropNulls col) of
         Left e -> throw e
         Right xs -> interQuartileRange' xs
 
--- | Calculates the Pearson's correlation coefficient between two given columns as a standalone value.
+{- | Calculates the Pearson's correlation coefficient between two given columns as a standalone value.
+Pairs with a null in either column are dropped.
+-}
 correlation :: T.Text -> T.Text -> DataFrame -> Maybe Double
 correlation first second df = do
-    f <- _getColumnAsDouble first df
-    s <- _getColumnAsDouble second df
+    let df' = filterJust first (filterJust second df)
+    f <- _getColumnAsDouble first df'
+    s <- _getColumnAsDouble second df'
     correlation' f s
 
 _getColumnAsDouble :: T.Text -> DataFrame -> Maybe (VU.Vector Double)
 _getColumnAsDouble name df = case getColumn name df of
-    Just (UnboxedColumn _ (f :: VU.Vector a)) -> case testEquality (typeRep @a) (typeRep @Double) of
-        Just Refl -> Just f
-        Nothing -> case sIntegral @a of
-            STrue -> Just (VU.map fromIntegral f)
-            SFalse -> case sFloating @a of
-                STrue -> Just (VU.map realToFrac f)
-                SFalse -> Nothing
+    Just col -> case dropNulls col of
+        UnboxedColumn _ (f :: VU.Vector a) ->
+            case testEquality (typeRep @a) (typeRep @Double) of
+                Just Refl -> Just f
+                Nothing -> case sIntegral @a of
+                    STrue -> Just (VU.map fromIntegral f)
+                    SFalse -> case sFloating @a of
+                        STrue -> Just (VU.map realToFrac f)
+                        SFalse -> Nothing
+        _ -> Nothing
     Nothing ->
         throw $
             ColumnsNotFoundException [name] "_getColumnAsDouble" (M.keys $ columnIndices df)
-    _ -> Nothing
 {-# INLINE _getColumnAsDouble #-}
 
 optionalToDoubleVector :: (Real a) => V.Vector (Maybe a) -> VU.Vector Double
@@ -265,17 +291,18 @@ sum ::
     forall a. (Columnable a, Num a) => Expr a -> DataFrame -> a
 sum (Col name) df = case getColumn name df of
     Nothing -> throw $ ColumnsNotFoundException [name] "sum" (M.keys $ columnIndices df)
-    Just ((UnboxedColumn _ (column :: VU.Vector a'))) -> case testEquality (typeRep @a') (typeRep @a) of
-        Just Refl -> VG.sum column
-        Nothing -> 0
-    Just ((BoxedColumn _ (column :: V.Vector a'))) -> case testEquality (typeRep @a') (typeRep @a) of
-        Just Refl -> VG.sum column
-        Nothing -> 0
-    Just (PackedText _ _) -> 0
-    Just (MergedColumn _ _) -> 0 -- matches the old eager These column (type never Num)
+    Just c -> case dropNulls c of
+        UnboxedColumn _ (column :: VU.Vector a') -> case testEquality (typeRep @a') (typeRep @a) of
+            Just Refl -> VG.sum column
+            Nothing -> 0
+        BoxedColumn _ (column :: V.Vector a') -> case testEquality (typeRep @a') (typeRep @a) of
+            Just Refl -> VG.sum column
+            Nothing -> 0
+        PackedText _ _ -> 0
+        MergedColumn _ _ -> 0 -- never numeric
 sum expr df = case interpret df expr of
     Left e -> throw e
-    Right (TColumn xs) -> case toVector @a @V.Vector xs of
+    Right (TColumn xs) -> case toVector @a @V.Vector (dropNulls xs) of
         Left e -> throw e
         Right xs' -> VG.sum xs'
 
@@ -416,7 +443,10 @@ summarize df =
 
 -- | Round a @Double@ to Specified Precision
 roundTo :: Int -> Double -> Double
-roundTo n x = fromInteger (round $ x * 10 ^ n) / 10.0 ^^ n
+roundTo n x
+    -- round on NaN is garbage
+    | isNaN x = x
+    | otherwise = fromInteger (round $ x * 10 ^ n) / 10.0 ^^ n
 
 toPct2dp :: Double -> String
 toPct2dp x
diff --git a/dataframe.cabal b/dataframe.cabal
index b7310c6f..5089bdf9 100644
--- a/dataframe.cabal
+++ b/dataframe.cabal
@@ -319,7 +319,6 @@ test-suite tests
                    Internal.DictEncode,
                    Internal.Markdown,
                    Internal.PackedText,
-                   Internal.Parsing,
                    PackedTextMigration,
                    PrettyPrint,
                    Learn.Denotation,
diff --git a/tests/Internal/Parsing.hs b/tests/Internal/Parsing.hs
deleted file mode 100644
index 4ad13e53..00000000
--- a/tests/Internal/Parsing.hs
+++ /dev/null
@@ -1,237 +0,0 @@
-{-# LANGUAGE OverloadedStrings #-}
-
-module Internal.Parsing where
-
-import DataFrame.Internal.Parsing
-import Test.HUnit
-
--- isNullish: recognized null strings
-
-isNullishEmptyString :: Test
-isNullishEmptyString =
-    TestCase (assertBool "empty string is nullish" (isNullish ""))
-
-isNullishNA :: Test
-isNullishNA = TestCase (assertBool "NA is nullish" (isNullish "NA"))
-
-isNullishNULL :: Test
-isNullishNULL = TestCase (assertBool "NULL is nullish" (isNullish "NULL"))
-
-isNullishNull :: Test
-isNullishNull = TestCase (assertBool "null is nullish" (isNullish "null"))
-
-isNullishNaN :: Test
-isNullishNaN = TestCase (assertBool "nan is nullish" (isNullish "nan"))
-
-isNullishNaNMixed :: Test
-isNullishNaNMixed = TestCase (assertBool "NaN is nullish" (isNullish "NaN"))
-
-isNullishNANUpper :: Test
-isNullishNANUpper = TestCase (assertBool "NAN is nullish" (isNullish "NAN"))
-
-isNullishNothing :: Test
-isNullishNothing =
-    TestCase (assertBool "Nothing is nullish" (isNullish "Nothing"))
-
-isNullishSpace :: Test
-isNullishSpace =
-    TestCase (assertBool "single space is nullish" (isNullish " "))
-
-isNullishNSlashA :: Test
-isNullishNSlashA = TestCase (assertBool "N/A is nullish" (isNullish "N/A"))
-
--- isNullish: values that are NOT null
-
-notNullishNumber :: Test
-notNullishNumber =
-    TestCase (assertBool "\"42\" is not nullish" (not (isNullish "42")))
-
-notNullishText :: Test
-notNullishText =
-    TestCase (assertBool "\"hello\" is not nullish" (not (isNullish "hello")))
-
-notNullishTrue :: Test
-notNullishTrue =
-    TestCase (assertBool "\"True\" is not nullish" (not (isNullish "True")))
-
-notNullishDouble :: Test
-notNullishDouble =
-    TestCase (assertBool "\"3.14\" is not nullish" (not (isNullish "3.14")))
-
-notNullishZero :: Test
-notNullishZero =
-    TestCase (assertBool "\"0\" is not nullish" (not (isNullish "0")))
-
--- readBool: positive cases
-
-readBoolTrue :: Test
-readBoolTrue =
-    TestCase (assertEqual "readBool \"True\"" (Just True) (readBool "True"))
-
-readBoolTrueLower :: Test
-readBoolTrueLower =
-    TestCase (assertEqual "readBool \"true\"" (Just True) (readBool "true"))
-
-readBoolTrueUpper :: Test
-readBoolTrueUpper =
-    TestCase (assertEqual "readBool \"TRUE\"" (Just True) (readBool "TRUE"))
-
-readBoolFalse :: Test
-readBoolFalse =
-    TestCase (assertEqual "readBool \"False\"" (Just False) (readBool "False"))
-
-readBoolFalseLower :: Test
-readBoolFalseLower =
-    TestCase (assertEqual "readBool \"false\"" (Just False) (readBool "false"))
-
-readBoolFalseUpper :: Test
-readBoolFalseUpper =
-    TestCase (assertEqual "readBool \"FALSE\"" (Just False) (readBool "FALSE"))
-
--- readBool: values that are not booleans
-
-readBoolDigit :: Test
-readBoolDigit =
-    TestCase (assertEqual "readBool \"1\" is Nothing" Nothing (readBool "1"))
-
-readBoolYes :: Test
-readBoolYes =
-    TestCase (assertEqual "readBool \"yes\" is Nothing" Nothing (readBool "yes"))
-
-readBoolEmpty :: Test
-readBoolEmpty =
-    TestCase (assertEqual "readBool \"\" is Nothing" Nothing (readBool ""))
-
-readBoolPartialTrue :: Test
-readBoolPartialTrue =
-    TestCase (assertEqual "readBool \"Tru\" is Nothing" Nothing (readBool "Tru"))
-
--- readInt
-
-readIntPositive :: Test
-readIntPositive =
-    TestCase (assertEqual "readInt \"42\"" (Just 42) (readInt "42"))
-
-readIntNegative :: Test
-readIntNegative =
-    TestCase (assertEqual "readInt \"-17\"" (Just (-17)) (readInt "-17"))
-
-readIntZero :: Test
-readIntZero =
-    TestCase (assertEqual "readInt \"0\"" (Just 0) (readInt "0"))
-
--- readInt strips whitespace before parsing
-readIntLeadingSpace :: Test
-readIntLeadingSpace =
-    TestCase
-        ( assertEqual
-            "readInt \" 5 \" (strips whitespace)"
-            (Just 5)
-            (readInt " 5 ")
-        )
-
-readIntFloat :: Test
-readIntFloat =
-    TestCase
-        (assertEqual "readInt \"3.14\" is Nothing" Nothing (readInt "3.14"))
-
-readIntText :: Test
-readIntText =
-    TestCase (assertEqual "readInt \"abc\" is Nothing" Nothing (readInt "abc"))
-
-readIntEmpty :: Test
-readIntEmpty =
-    TestCase (assertEqual "readInt \"\" is Nothing" Nothing (readInt ""))
-
--- trailing non-digits must make the parse fail
-readIntPartialSuffix :: Test
-readIntPartialSuffix =
-    TestCase
-        ( assertEqual
-            "readInt \"42abc\" is Nothing"
-            Nothing
-            (readInt "42abc")
-        )
-
--- readDouble
-
-readDoublePositive :: Test
-readDoublePositive =
-    TestCase
-        (assertEqual "readDouble \"3.14\"" (Just 3.14) (readDouble "3.14"))
-
-readDoubleNegative :: Test
-readDoubleNegative =
-    TestCase
-        (assertEqual "readDouble \"-1.5\"" (Just (-1.5)) (readDouble "-1.5"))
-
-readDoubleWholeNumber :: Test
-readDoubleWholeNumber =
-    TestCase
-        ( assertEqual
-            "readDouble \"42\" parses as 42.0"
-            (Just 42.0)
-            (readDouble "42")
-        )
-
-readDoubleText :: Test
-readDoubleText =
-    TestCase
-        (assertEqual "readDouble \"abc\" is Nothing" Nothing (readDouble "abc"))
-
-readDoubleEmpty :: Test
-readDoubleEmpty =
-    TestCase
-        (assertEqual "readDouble \"\" is Nothing" Nothing (readDouble ""))
-
-readDoublePartialSuffix :: Test
-readDoublePartialSuffix =
-    TestCase
-        ( assertEqual
-            "readDouble \"3.14abc\" is Nothing"
-            Nothing
-            (readDouble "3.14abc")
-        )
-
-tests :: [Test]
-tests =
-    [ TestLabel "isNullishEmptyString" isNullishEmptyString
-    , TestLabel "isNullishNA" isNullishNA
-    , TestLabel "isNullishNULL" isNullishNULL
-    , TestLabel "isNullishNull" isNullishNull
-    , TestLabel "isNullishNaN" isNullishNaN
-    , TestLabel "isNullishNaNMixed" isNullishNaNMixed
-    , TestLabel "isNullishNANUpper" isNullishNANUpper
-    , TestLabel "isNullishNothing" isNullishNothing
-    , TestLabel "isNullishSpace" isNullishSpace
-    , TestLabel "isNullishNSlashA" isNullishNSlashA
-    , TestLabel "notNullishNumber" notNullishNumber
-    , TestLabel "notNullishText" notNullishText
-    , TestLabel "notNullishTrue" notNullishTrue
-    , TestLabel "notNullishDouble" notNullishDouble
-    , TestLabel "notNullishZero" notNullishZero
-    , TestLabel "readBoolTrue" readBoolTrue
-    , TestLabel "readBoolTrueLower" readBoolTrueLower
-    , TestLabel "readBoolTrueUpper" readBoolTrueUpper
-    , TestLabel "readBoolFalse" readBoolFalse
-    , TestLabel "readBoolFalseLower" readBoolFalseLower
-    , TestLabel "readBoolFalseUpper" readBoolFalseUpper
-    , TestLabel "readBoolDigit" readBoolDigit
-    , TestLabel "readBoolYes" readBoolYes
-    , TestLabel "readBoolEmpty" readBoolEmpty
-    , TestLabel "readBoolPartialTrue" readBoolPartialTrue
-    , TestLabel "readIntPositive" readIntPositive
-    , TestLabel "readIntNegative" readIntNegative
-    , TestLabel "readIntZero" readIntZero
-    , TestLabel "readIntLeadingSpace" readIntLeadingSpace
-    , TestLabel "readIntFloat" readIntFloat
-    , TestLabel "readIntText" readIntText
-    , TestLabel "readIntEmpty" readIntEmpty
-    , TestLabel "readIntPartialSuffix" readIntPartialSuffix
-    , TestLabel "readDoublePositive" readDoublePositive
-    , TestLabel "readDoubleNegative" readDoubleNegative
-    , TestLabel "readDoubleWholeNumber" readDoubleWholeNumber
-    , TestLabel "readDoubleText" readDoubleText
-    , TestLabel "readDoubleEmpty" readDoubleEmpty
-    , TestLabel "readDoublePartialSuffix" readDoublePartialSuffix
-    ]
diff --git a/tests/Learn/MetricsTests.hs b/tests/Learn/MetricsTests.hs
index d43ccea6..87089b05 100644
--- a/tests/Learn/MetricsTests.hs
+++ b/tests/Learn/MetricsTests.hs
@@ -1,8 +1,11 @@
 {-# LANGUAGE OverloadedStrings #-}
+{-# LANGUAGE ScopedTypeVariables #-}
 {-# LANGUAGE TypeApplications #-}
 
 module Learn.MetricsTests (tests) where
 
+import qualified Control.Exception as E
+
 import qualified DataFrame as D
 import qualified DataFrame.Functions as F
 import qualified DataFrame.Internal.Column as DI
@@ -43,6 +46,16 @@ testRegressionMetrics = TestCase $ do
     assertBool "rmse" (close 1e-9 (rmse p t) 0.5)
     assertBool "mae" (close 1e-9 (mae p t) 0.25)
     assertBool "r2 in range" (r2 p t <= 1)
+    assertBool
+        "mse averages over compared pairs"
+        (close 1e-9 (mse (VU.fromList [2, 2]) (VU.fromList [0, 0, 0, 0])) 4)
+    assertBool
+        "mae averages over compared pairs"
+        (close 1e-9 (mae (VU.fromList [2, 2]) (VU.fromList [0, 0, 0, 0])) 2)
+    r <- E.try (E.evaluate (mse VU.empty (VU.fromList [5, 5, 5])))
+    case r of
+        Left (_ :: E.SomeException) -> pure ()
+        Right v -> assertFailure ("mse with no pairs returned " ++ show v)
 
 testMulticlassMetrics :: Test
 testMulticlassMetrics = TestCase $ do
diff --git a/tests/Main.hs b/tests/Main.hs
index 7689e2f7..3263a53d 100644
--- a/tests/Main.hs
+++ b/tests/Main.hs
@@ -17,7 +17,6 @@ import qualified Internal.ColumnBuilder
 import qualified Internal.DictEncode
 import qualified Internal.Markdown
 import qualified Internal.PackedText
-import qualified Internal.Parsing
 import qualified LazyParity
 import qualified LazyParquet
 import qualified LazyProjection
@@ -75,7 +74,6 @@ tests =
             ++ Internal.DictEncode.tests
             ++ Internal.Markdown.tests
             ++ Internal.PackedText.tests
-            ++ Internal.Parsing.tests
             ++ Learn.Denotation.tests
             ++ Learn.Models.tests
             ++ Learn.TypedModel.tests
diff --git a/tests/Operations/ParallelGroupBy.hs b/tests/Operations/ParallelGroupBy.hs
index 4b5b9310..c808387b 100644
--- a/tests/Operations/ParallelGroupBy.hs
+++ b/tests/Operations/ParallelGroupBy.hs
@@ -109,7 +109,11 @@ aggParityFor n =
                 ]
             seqDf = D.aggregate aggs (groupBySeq ["ki", "kt"] df)
             parDf = D.aggregate aggs (groupByPar ["ki", "kt"] df)
-         in assertEqual ("aggregate parity n=" ++ show n) seqDf parDf
+         in -- render: NaN /= NaN under Eq
+            assertEqual
+                ("aggregate parity n=" ++ show n)
+                (D.toMarkdown seqDf)
+                (D.toMarkdown parDf)
 
 collisionParity :: Test
 collisionParity =
diff --git a/tests/Operations/Statistics.hs b/tests/Operations/Statistics.hs
index db907654..ae253d8f 100644
--- a/tests/Operations/Statistics.hs
+++ b/tests/Operations/Statistics.hs
@@ -6,6 +6,7 @@ module Operations.Statistics where
 
 import qualified Data.Vector.Unboxed as VU
 import qualified DataFrame as D
+import qualified DataFrame.Functions as F
 import qualified DataFrame.Internal.Column as DI
 import qualified DataFrame.Internal.Statistics as D
 
@@ -56,6 +57,7 @@ skewnessOfSymmetricDataSet =
             0
         )
 
+-- g1, matching scipy.stats.skew
 skewnessOfSimpleDataSet :: Test
 skewnessOfSimpleDataSet =
     TestCase
@@ -63,7 +65,7 @@ skewnessOfSimpleDataSet =
             "Skewness of a simple data set"
             ( abs
                 ( D.skewness' (VU.fromList [25 :: Int, 28, 26, 30, 40, 50, 40])
-                    - 0.566_731_633_676
+                    - 0.612_140_127_240_396_6
                 )
                 < 1e-12
             )
@@ -254,6 +256,167 @@ correlationMissingColumn =
             (print $ D.correlation "x" "missingcol" correlationDf)
         )
 
+-- stats must skip null slots
+
+nullableDf :: D.DataFrame
+nullableDf =
+    D.fromNamedColumns [("x", DI.fromList [Just (10 :: Double), Nothing, Just 20])]
+
+nullableIntDf :: D.DataFrame
+nullableIntDf =
+    D.fromNamedColumns [("n", DI.fromList [Just (10 :: Int), Nothing, Just 20])]
+
+nullableBoxedDf :: D.DataFrame
+nullableBoxedDf =
+    D.fromNamedColumns [("b", DI.fromList [Just (10 :: Integer), Nothing, Just 20])]
+
+meanIgnoresNulls :: Test
+meanIgnoresNulls =
+    TestCase
+        (assertEqual "mean skips nulls" 15.0 (D.mean (F.col @Double "x") nullableDf))
+
+meanExprIgnoresNulls :: Test
+meanExprIgnoresNulls =
+    TestCase
+        ( assertEqual
+            "mean over a derived nullable expression skips nulls"
+            30.0
+            (D.mean (F.lift (* 2) (F.col @Double "x")) nullableDf)
+        )
+
+medianIgnoresNulls :: Test
+medianIgnoresNulls =
+    TestCase
+        (assertEqual "median skips nulls" 15.0 (D.median (F.col @Double "x") nullableDf))
+
+percentileIgnoresNulls :: Test
+percentileIgnoresNulls =
+    TestCase
+        ( assertEqual
+            "percentile skips nulls"
+            15.0
+            (D.percentile 50 (F.col @Double "x") nullableDf)
+        )
+
+stdDevIgnoresNulls :: Test
+stdDevIgnoresNulls =
+    TestCase
+        ( assertBool
+            "standard deviation skips nulls"
+            ( abs (D.standardDeviation (F.col @Double "x") nullableDf - 7.0710678118654755)
+                < 1e-12
+            )
+        )
+
+varianceIgnoresNulls :: Test
+varianceIgnoresNulls =
+    TestCase
+        ( assertEqual
+            "variance skips nulls"
+            50.0
+            (D.variance (F.col @Double "x") nullableDf)
+        )
+
+varianceExprIgnoresNulls :: Test
+varianceExprIgnoresNulls =
+    TestCase
+        ( assertEqual
+            "variance over a derived nullable expression skips nulls"
+            200.0
+            (D.variance (F.lift (* 2) (F.col @Double "x")) nullableDf)
+        )
+
+iqrIgnoresNulls :: Test
+iqrIgnoresNulls =
+    TestCase
+        ( assertEqual
+            "inter-quartile range skips nulls"
+            5.0
+            (D.interQuartileRange (F.col @Double "x") nullableDf)
+        )
+
+skewnessIgnoresNulls :: Test
+skewnessIgnoresNulls =
+    TestCase
+        ( let skewDf =
+                D.fromNamedColumns
+                    [("s", DI.fromList [Just (10 :: Double), Nothing, Just 20, Just 100, Just 11])]
+           in assertBool
+                "skewness skips nulls"
+                ( abs
+                    ( D.skewness (F.col @Double "s") skewDf
+                        - D.skewness' (VU.fromList [10 :: Double, 20, 100, 11])
+                    )
+                    < 1e-12
+                )
+        )
+
+genericPercentileIgnoresNulls :: Test
+genericPercentileIgnoresNulls =
+    TestCase
+        ( assertEqual
+            "genericPercentile skips the sentinel"
+            10
+            (D.genericPercentile 10 (F.col @Int "n") nullableIntDf)
+        )
+
+genericPercentileBoxedNullableDoesNotThrow :: Test
+genericPercentileBoxedNullableDoesNotThrow =
+    TestCase
+        ( assertEqual
+            "genericPercentile on a boxed nullable column skips the error thunk"
+            20
+            (D.genericPercentile 100 (F.col @Integer "b") nullableBoxedDf)
+        )
+
+genericPercentileMaybeViewKeepsNothing :: Test
+genericPercentileMaybeViewKeepsNothing =
+    TestCase
+        ( assertEqual
+            "a Maybe-typed view still sees its Nothings"
+            (Nothing :: Maybe Int)
+            (D.genericPercentile 0 (F.col @(Maybe Int) "n") nullableIntDf)
+        )
+
+sumUnboxedNullable :: Test
+sumUnboxedNullable =
+    TestCase
+        (assertEqual "sum skips null slots" 30.0 (D.sum (F.col @Double "x") nullableDf))
+
+sumBoxedNullableDoesNotThrow :: Test
+sumBoxedNullableDoesNotThrow =
+    TestCase
+        ( assertEqual
+            "sum on a boxed nullable column skips the error thunk"
+            (30 :: Integer)
+            (D.sum (F.col @Integer "b") nullableBoxedDf)
+        )
+
+correlationIgnoresNullRows :: Test
+correlationIgnoresNullRows =
+    TestCase
+        ( let dfc =
+                D.fromNamedColumns
+                    [ ("a", DI.fromList [Just (1 :: Double), Nothing, Just 3])
+                    , ("c", DI.fromList [1 :: Double, 2, 3])
+                    ]
+           in case D.correlation "a" "c" dfc of
+                Nothing -> assertFailure "Expected Just 1.0, got Nothing"
+                Just r ->
+                    assertBool
+                        "null rows are dropped pairwise"
+                        (abs (r - 1.0) < 1e-10)
+        )
+
+frequenciesNullableAsMaybe :: Test
+frequenciesNullableAsMaybe =
+    TestCase
+        ( assertEqual
+            "nulls are a Nothing category"
+            4 -- Statistic, Nothing, Just 10, Just 20
+            (D.nColumns (D.frequencies (F.col @(Maybe Int) "n") nullableIntDf))
+        )
+
 tests :: [Test]
 tests =
     [ TestLabel "medianOfOddLengthDataSet" medianOfOddLengthDataSet
@@ -281,4 +444,24 @@ tests =
     , TestLabel "correlationPerfectNegative" correlationPerfectNegative
     , TestLabel "correlationSelfIdentity" correlationSelfIdentity
     , TestLabel "correlationMissingColumn" correlationMissingColumn
+    , TestLabel "meanIgnoresNulls" meanIgnoresNulls
+    , TestLabel "meanExprIgnoresNulls" meanExprIgnoresNulls
+    , TestLabel "medianIgnoresNulls" medianIgnoresNulls
+    , TestLabel "percentileIgnoresNulls" percentileIgnoresNulls
+    , TestLabel "stdDevIgnoresNulls" stdDevIgnoresNulls
+    , TestLabel "varianceIgnoresNulls" varianceIgnoresNulls
+    , TestLabel "varianceExprIgnoresNulls" varianceExprIgnoresNulls
+    , TestLabel "iqrIgnoresNulls" iqrIgnoresNulls
+    , TestLabel "skewnessIgnoresNulls" skewnessIgnoresNulls
+    , TestLabel "genericPercentileIgnoresNulls" genericPercentileIgnoresNulls
+    , TestLabel
+        "genericPercentileBoxedNullableDoesNotThrow"
+        genericPercentileBoxedNullableDoesNotThrow
+    , TestLabel
+        "genericPercentileMaybeViewKeepsNothing"
+        genericPercentileMaybeViewKeepsNothing
+    , TestLabel "sumUnboxedNullable" sumUnboxedNullable
+    , TestLabel "sumBoxedNullableDoesNotThrow" sumBoxedNullableDoesNotThrow
+    , TestLabel "correlationIgnoresNullRows" correlationIgnoresNullRows
+    , TestLabel "frequenciesNullableAsMaybe" frequenciesNullableAsMaybe
     ]
diff --git a/tests/Operations/VectorKernel.hs b/tests/Operations/VectorKernel.hs
index b9c898cb..802e015d 100644
--- a/tests/Operations/VectorKernel.hs
+++ b/tests/Operations/VectorKernel.hs
@@ -123,7 +123,11 @@ parityCase n keys aggs =
             fast = D.aggregate aggs gdf
             ref = interpretOnly aggs gdf
             label = "n=" ++ show n ++ " keys=" ++ show keys ++ " #aggs=" ++ show (length aggs)
-         in assertEqual ("kernel==interpreter " ++ label) ref fast
+         in -- render: NaN /= NaN under Eq
+            assertEqual
+                ("kernel==interpreter " ++ label)
+                (D.toMarkdown ref)
+                (D.toMarkdown fast)
 
 tests :: [Test]
 tests =