diff --git a/mealy.cabal b/mealy.cabal
--- a/mealy.cabal
+++ b/mealy.cabal
@@ -1,6 +1,6 @@
 cabal-version: 2.4
 name:          mealy
-version:       0.0.1
+version:       0.0.2
 synopsis: See readme.md
 description: See readme.md for description.
 category: project
@@ -27,10 +27,11 @@
     generic-lens >= 2.0,
     hmatrix >= 0.20,
     lens,
+    matrix >= 0.3.6 && < 0.3.7,
     mwc-probability,
     mwc-random,
-    numhask >= 0.6 && < 0.7,
-    numhask-array >= 0.7 && < 0.8,
+    numhask >= 0.7 && < 0.8,
+    numhask-array >= 0.8 && < 0.9,
     primitive >= 0.7,
     profunctors >= 5.5,
     tdigest,
@@ -54,6 +55,8 @@
     -Wincomplete-record-updates
     -Wincomplete-uni-patterns
     -Wredundant-constraints
+    -fwrite-ide-info
+    -hiedir=.hie
 
 test-suite test
   type: exitcode-stdio-1.0
@@ -63,7 +66,7 @@
   build-depends:
     base >=4.7 && <5,
     doctest,
-    numhask >= 0.6 && < 0.7,
+    numhask >= 0.7 && < 0.8,
     mealy
   default-language: Haskell2010
   default-extensions:
@@ -77,3 +80,5 @@
     -Wincomplete-record-updates
     -Wincomplete-uni-patterns
     -Wredundant-constraints
+    -fwrite-ide-info
+    -hiedir=.hie
diff --git a/src/Data/Mealy.hs b/src/Data/Mealy.hs
--- a/src/Data/Mealy.hs
+++ b/src/Data/Mealy.hs
@@ -13,14 +13,15 @@
 {-# LANGUAGE PatternSynonyms #-}
 {-# LANGUAGE RankNTypes #-}
 {-# LANGUAGE ScopedTypeVariables #-}
+{-# LANGUAGE StrictData #-}
 {-# LANGUAGE TupleSections #-}
-{-# LANGUAGE TypeOperators #-}
 {-# LANGUAGE NoImplicitPrelude #-}
 {-# OPTIONS_GHC -Wall #-}
+{-# OPTIONS_GHC -Wno-incomplete-patterns #-}
+{-# OPTIONS_GHC -Wno-incomplete-uni-patterns #-}
 {-# OPTIONS_GHC -Wno-name-shadowing #-}
 {-# OPTIONS_GHC -Wno-orphans #-}
 {-# OPTIONS_GHC -Wno-type-defaults #-}
-{-# OPTIONS_GHC -Wno-incomplete-patterns #-}
 
 -- | Online statistics for ordered data (such as time-series data), modelled as [mealy machines](https://en.wikipedia.org/wiki/Mealy_machine)
 module Data.Mealy
@@ -72,18 +73,21 @@
   )
 where
 
-import Control.Lens hiding ((:>), Empty, Unwrapped, Wrapped, index, (|>))
+import Control.Lens hiding (Empty, Unwrapped, Wrapped, index, (:>), (|>))
 import Data.Fold hiding (M)
 import Data.Functor.Rep
 import Data.Generics.Labels ()
+-- import qualified NumHask.Array.HMatrix as HM
+
+-- import qualified Numeric.LinearAlgebra as LA
+
+import qualified Data.Matrix as M
 import qualified Data.Sequence as Seq
 import qualified NumHask.Array.Fixed as F
-import qualified NumHask.Array.HMatrix as HM
 import NumHask.Array.Shape (HasShape)
 import NumHask.Prelude hiding (L1, State, StateT, asum, fold, get, replace, runState, runStateT, state)
-import qualified Numeric.Backprop as B
 import Numeric.Backprop (BVar, Reifies, W)
-import qualified Numeric.LinearAlgebra as LA
+import qualified Numeric.Backprop as B
 import qualified Prelude.Backprop as PB
 import qualified Prelude as P
 
@@ -104,23 +108,21 @@
 -- >>> xs2 <- rvs g 10000
 -- >>> xsp <- rvsp g 10000 0.8
 
-{- | A 'Mealy' is a triple of functions
-
- * (a -> b) __inject__ Convert an input into the state type.
- * (b -> a -> b) __step__ Update state given prior state and (new) input.
- * (c -> b) __extract__ Convert state to the output type.
-
- By adopting this order, a Mealy sum looks like:
-
-> M id (+) id
-
-where the first id is the initial injection to a contravariant position, and the second id is the covriant extraction.
-
- __inject__ kicks off state on the initial element of the Foldable, but is otherwise be independent of __step__.
-
-> scan (M e s i) (x : xs) = e <$> scanl' s (i x) xs
-
--}
+-- | A 'Mealy' is a triple of functions
+--
+-- * (a -> b) __inject__ Convert an input into the state type.
+-- * (b -> a -> b) __step__ Update state given prior state and (new) input.
+-- * (c -> b) __extract__ Convert state to the output type.
+--
+-- By adopting this order, a Mealy sum looks like:
+--
+-- > M id (+) id
+--
+-- where the first id is the initial injection to a contravariant position, and the second id is the covriant extraction.
+--
+-- __inject__ kicks off state on the initial element of the Foldable, but is otherwise be independent of __step__.
+--
+-- > scan (M e s i) (x : xs) = e <$> scanl' s (i x) xs
 newtype Mealy a b = Mealy {l1 :: L1 a b}
   deriving (Profunctor, Category) via L1
   deriving (Functor, Applicative) via L1 a
@@ -148,10 +150,9 @@
 scan (M i s e) (x : xs) = fromList (e <$> scanl' s (i x) xs)
 
 -- | Most common statistics are averages, which are some sort of aggregation of values (sum) and some sort of sample size (count).
-newtype Averager a b
-  = Averager
-      { sumCount :: (a, b)
-      }
+newtype Averager a b = Averager
+  { sumCount :: (a, b)
+  }
   deriving (Eq, Show)
 
 -- | Pattern for an 'Averager'.
@@ -216,7 +217,7 @@
 --
 -- >>> fold (absma 1) xs0
 -- 0.7894201075535578
-absma :: (Divisive a, Additive a, Signed a) => a -> Mealy a a
+absma :: (Divisive a, Signed a) => a -> Mealy a a
 absma r = online abs (* r)
 {-# INLINEABLE absma #-}
 
@@ -330,13 +331,12 @@
 reg1 :: (ExpField a) => Mealy a a -> Mealy (a, a) (a, a)
 reg1 m = (,) <$> alpha1 m <*> beta1 m
 
-data RegressionState (n :: Nat) a
-  = RegressionState
-      { _xx :: F.Array '[n, n] a,
-        _x :: F.Array '[n] a,
-        _xy :: F.Array '[n] a,
-        _y :: a
-      }
+data RegressionState (n :: Nat) a = RegressionState
+  { _xx :: F.Array '[n, n] a,
+    _x :: F.Array '[n] a,
+    _xy :: F.Array '[n] a,
+    _y :: a
+  }
   deriving (Functor)
 
 -- | multiple regression
@@ -357,22 +357,37 @@
 -- > let zs = zip (zipWith (\x y -> fromList [x,y] :: F.Array '[2] Double) xs1 xs2) ys
 -- > fold (beta 0.99) zs
 -- [0.4982692361226971, 1.038192474255091]
-beta :: (Field a, LA.Field a, KnownNat n) => a -> Mealy (F.Array '[n] a, a) (F.Array '[n] a)
+beta :: (Field a, KnownNat n, Fractional a, Eq a) => a -> Mealy (F.Array '[n] a, a) (F.Array '[n] a)
 beta r = M inject step extract
   where
+    -- extract :: Averager (RegressionState n a) a -> (F.Array '[n] a)
     extract (A (RegressionState xx x xy y) c) =
-      liftHM2
-        (\a b -> LA.pinv a LA.<> LA.tr b)
-        ((one / c) *. (xx - F.expand (*) x x))
-        ((xy - (y *. x)) .* (one / c))
+      (\a b -> inverse a `F.mult` b)
+        ((one / c) .* (xx - F.expand (*) x x))
+        ((xy - (y .* x)) *. (one / c))
     step x (xs, y) = rsOnline r x (inject (xs, y))
+    -- inject :: (F.Array '[n] a, a) -> Averager (RegressionState n a) a
     inject (xs, y) =
-      A (RegressionState (F.expand (*) xs xs) xs (y *. xs) y) one
+      A (RegressionState (F.expand (*) xs xs) xs (y .* xs) y) one
 {-# INLINEABLE beta #-}
 
-liftHM2 :: (LA.Field a, HasShape s, HasShape s', HasShape s'') => (LA.Matrix a -> LA.Matrix a -> LA.Matrix a) -> F.Array s a -> F.Array s' a -> F.Array s'' a
-liftHM2 f a b = HM.toFixed $ HM.Array $ f (HM.unArray . HM.fromFixed $ a) (HM.unArray . HM.fromFixed $ b)
+toMatrix :: (KnownNat n, KnownNat m) => F.Array [m, n] a -> M.Matrix a
+toMatrix a = M.matrix m n (index a . (\(i, j) -> [i, j]))
+  where
+    (m : n : _) = F.shape a
 
+fromMatrix :: (KnownNat n, KnownNat m) => M.Matrix a -> F.Array [m, n] a
+fromMatrix = fromList . M.toList
+
+data MatrixException = MatrixException
+  deriving (Show)
+
+instance Exception MatrixException
+
+-- | The inverse of a square matrix.
+inverse :: (KnownNat n, Fractional a, Eq a) => F.Array [n, n] a -> F.Array [n, n] a
+inverse = either (const $ throw MatrixException) fromMatrix . M.inverse . toMatrix
+
 rsOnline :: (Field a, KnownNat n) => a -> Averager (RegressionState n a) a -> Averager (RegressionState n a) a -> Averager (RegressionState n a) a
 rsOnline r (A (RegressionState xx x xy y) c) (A (RegressionState xx' x' xy' y') c') =
   A (RegressionState (liftR2 d xx xx') (liftR2 d x x') (liftR2 d xy xy') (d y y')) (d c c')
@@ -380,7 +395,7 @@
     d s s' = r * s + s'
 
 -- | alpha in a multiple regression
-alpha :: (LA.Field a, ExpField a, KnownNat n) => a -> Mealy (F.Array '[n] a, a) a
+alpha :: (ExpField a, KnownNat n, Fractional a, Eq a) => a -> Mealy (F.Array '[n] a, a) a
 alpha r = (\xs b y -> y - sum (liftR2 (*) b xs)) <$> lmap fst (arrayify $ ma r) <*> beta r <*> lmap snd (ma r)
 {-# INLINEABLE alpha #-}
 
@@ -397,7 +412,7 @@
 -- > let zs = zip (zipWith (\x y -> fromList [x,y] :: F.Array '[2] Double) xs1 xs2) ys
 -- > fold (reg 0.99) zs
 -- ([0.4982692361226971, 1.038192474255091],2.087160803386695e-3)
-reg :: (LA.Field a, ExpField a, KnownNat n) => a -> Mealy (F.Array '[n] a, a) (F.Array '[n] a, a)
+reg :: (ExpField a, KnownNat n, Fractional a, Eq a) => a -> Mealy (F.Array '[n] a, a) (F.Array '[n] a, a)
 reg r = (,) <$> beta r <*> alpha r
 {-# INLINEABLE reg #-}
 
@@ -497,15 +512,14 @@
 -- s_{t+1} & = (alpha_t^s + beta_t^{x->s} * ma_t^x + beta_t^{s->s} * std_t^x) * N(0,1)
 -- \end{align}
 -- \]
-data Model1
-  = Model1
-      { alphaX :: Double,
-        alphaS :: Double,
-        betaMa2X :: Double,
-        betaMa2S :: Double,
-        betaStd2X :: Double,
-        betaStd2S :: Double
-      }
+data Model1 = Model1
+  { alphaX :: Double,
+    alphaS :: Double,
+    betaMa2X :: Double,
+    betaMa2S :: Double,
+    betaStd2X :: Double,
+    betaStd2S :: Double
+  }
   deriving (Eq, Show, Generic)
 
 zeroModel1 :: Model1
@@ -528,8 +542,10 @@
               + (m1 ^. #betaStd2S) * (s - 1)
           )
         + m1 ^. #alphaX
-        + (m1 ^. #betaMa2X) * m
-        + (m1 ^. #betaStd2X) * (s - 1)
+        + (m1 ^. #betaMa2X)
+        * m
+        + (m1 ^. #betaStd2X)
+        * (s - 1)
 
 foldB :: (Reifies s W) => (BVar s Double -> BVar s Double) -> BVar s Double -> BVar s [Double] -> BVar s Double
 foldB f r xs = divide (PB.foldl' (step' f r) (B.T2 0 0) xs)
@@ -542,12 +558,11 @@
 
 -- | A rough Median.
 -- The average absolute value of the stat is used to callibrate estimate drift towards the median
-data Medianer a b
-  = Medianer
-      { medAbsSum :: a,
-        medCount :: b,
-        medianEst :: a
-      }
+data Medianer a b = Medianer
+  { medAbsSum :: a,
+    medCount :: b,
+    medianEst :: a
+  }
 
 -- | onlineL1' takes a function and turns it into a `Mealy` where the step is an incremental update of an (isomorphic) median statistic.
 onlineL1' ::
diff --git a/src/Data/Mealy/Quantiles.hs b/src/Data/Mealy/Quantiles.hs
--- a/src/Data/Mealy/Quantiles.hs
+++ b/src/Data/Mealy/Quantiles.hs
@@ -1,5 +1,6 @@
-{-# LANGUAGE NoImplicitPrelude #-}
 {-# LANGUAGE DataKinds #-}
+{-# LANGUAGE RebindableSyntax #-}
+{-# LANGUAGE StrictData #-}
 
 module Data.Mealy.Quantiles
   ( median,
@@ -17,12 +18,11 @@
 import qualified Data.Vector.Unboxed as VU
 import NumHask.Prelude hiding (fold)
 
-data OnlineTDigest
-  = OnlineTDigest
-      { td :: TDigest 25,
-        tdN :: Int,
-        tdRate :: Double
-      }
+data OnlineTDigest = OnlineTDigest
+  { td :: TDigest 25,
+    tdN :: Int,
+    tdRate :: Double
+  }
   deriving (Show)
 
 emptyOnlineTDigest :: Double -> OnlineTDigest