diff --git a/ChangeLog.md b/ChangeLog.md
--- a/ChangeLog.md
+++ b/ChangeLog.md
@@ -1,5 +1,13 @@
 # Revision history for rhine-gloss
 
+
+## 1.3
+
+* Dropped `dunai` dependency in favour of state automata.
+  See [the versions readme](./versions.md) for details.
+* Support GHC 9.6 and 9.8
+* Updated to monad-bayes 1.3
+
 ## 1.2
 
 * Updated to monad-bayes 1.2.0
diff --git a/app/Main.hs b/app/Main.hs
--- a/app/Main.hs
+++ b/app/Main.hs
@@ -39,8 +39,8 @@
 import Control.Monad.Bayes.Population hiding (hoist)
 import Control.Monad.Bayes.Sampler.Strict
 
--- dunai
-import Control.Monad.Trans.MSF.Except
+-- automaton
+import Data.Automaton.Trans.Except
 
 -- rhine
 import FRP.Rhine
@@ -171,7 +171,7 @@
 
 -- | The number of particles used in the filter. Change according to available computing power.
 nParticles :: Int
-nParticles = 100
+nParticles = 400
 
 -- * Visualization
 
@@ -239,21 +239,11 @@
   arrMCl paintIO -< toThermometer $ translate 0 (double2Float temperature * thermometerScale) $ color (withAlpha (double2Float $ exp $ 0.2 * ln probability) white) $ rectangleSolid thermometerWidth 2
 
 drawParticles :: BehaviourF App td [(Pos, Log Double)] ()
-drawParticles = proc particlesPosition -> do
-  case particlesPosition of
-    [] -> returnA -< ()
-    p : ps -> do
-      drawParticle -< p
-      drawParticles -< ps
+drawParticles = traverseS_ drawParticle
 
 -- FIXME abstract using a library
 drawParticlesTemperature :: BehaviourF App td [(Temperature, Log Double)] ()
-drawParticlesTemperature = proc particlesPosition -> do
-  case particlesPosition of
-    [] -> returnA -< ()
-    p : ps -> do
-      drawParticleTemperature -< p
-      drawParticlesTemperature -< ps
+drawParticlesTemperature = traverseS_ drawParticleTemperature
 
 glossSettings :: GlossSettings
 glossSettings =
@@ -380,7 +370,7 @@
     }
 
 {- | The part of the program which simulates latent position and sensor,
-   running 100 times a second.
+   running 10 times a second.
 -}
 modelRhine :: Rhine (GlossConcT IO) (LiftClock IO GlossConcT (Millisecond 100)) Temperature (Temperature, (Sensor, Pos))
 modelRhine = hoistClSF sampleIOGloss (clId &&& genModelWithoutTemperature) @@ liftClock waitClock
@@ -398,19 +388,19 @@
 -}
 inference :: Rhine (GlossConcT IO) (LiftClock IO GlossConcT Busy) (Temperature, (Sensor, Pos)) Result
 inference = hoistClSF sampleIOGloss inferenceBehaviour @@ liftClock Busy
-  where
-    inferenceBehaviour :: (MonadDistribution m, Diff td ~ Double, MonadIO m) => BehaviourF m td (Temperature, (Sensor, Pos)) Result
-    inferenceBehaviour = proc (temperature, (measured, latent)) -> do
-      positionsAndTemperatures <- runPopulationCl nParticles resampleSystematic posteriorTemperatureProcess -< measured
-      returnA
-        -<
-          Result
-            { temperature
-            , measured
-            , latent
-            , particlesPosition = first snd <$> positionsAndTemperatures
-            , particlesTemperature = first fst <$> positionsAndTemperatures
-            }
+
+inferenceBehaviour :: (MonadDistribution m, Diff td ~ Double, MonadIO m) => BehaviourF m td (Temperature, (Sensor, Pos)) Result
+inferenceBehaviour = proc (temperature, (measured, latent)) -> do
+  positionsAndTemperatures <- runPopulationCl nParticles resampleSystematic posteriorTemperatureProcess -< measured
+  returnA
+    -<
+      Result
+        { temperature
+        , measured
+        , latent
+        , particlesPosition = first snd <$> positionsAndTemperatures
+        , particlesTemperature = first fst <$> positionsAndTemperatures
+        }
 
 -- | Visualize the current 'Result' at a rate controlled by the @gloss@ backend, usually 30 FPS.
 visualisationRhine :: Rhine (GlossConcT IO) (GlossClockUTC GlossSimClockIO) Result ()
diff --git a/rhine-bayes.cabal b/rhine-bayes.cabal
--- a/rhine-bayes.cabal
+++ b/rhine-bayes.cabal
@@ -1,42 +1,47 @@
-name:                rhine-bayes
-version:             1.2
-synopsis:            monad-bayes backend for Rhine
+name: rhine-bayes
+version: 1.3
+synopsis: monad-bayes backend for Rhine
 description:
   This package provides a backend to the @monad-bayes@ library,
   enabling you to write stochastic processes as signal functions,
   and performing online machine learning on them.
-license:             BSD3
-license-file:        LICENSE
-author:              Manuel Bärenz
-maintainer:          programming@manuelbaerenz.de
+
+license: BSD3
+license-file: LICENSE
+author: Manuel Bärenz
+maintainer: programming@manuelbaerenz.de
 -- copyright:
-category:            FRP
-build-type:          Simple
-extra-doc-files:     README.md ChangeLog.md
-cabal-version:       2.0
+category: FRP
+build-type: Simple
+extra-doc-files:
+  ChangeLog.md
+  README.md
 
+cabal-version: 2.0
+
 source-repository head
-  type:     git
+  type: git
   location: git@github.com:turion/rhine.git
 
 source-repository this
-  type:     git
+  type: git
   location: git@github.com:turion/rhine.git
-  tag:      v1.1
+  tag: v1.3
 
 library
-  exposed-modules:
-    FRP.Rhine.Bayes
-  other-modules:
-    Data.MonadicStreamFunction.Bayes
-  build-depends:       base         >= 4.11 && < 4.18
-                     , transformers >= 0.5
-                     , rhine        == 1.2
-                     , dunai        ^>= 0.11
-                     , log-domain   >= 0.12
-                     , monad-bayes  ^>= 1.2
-  hs-source-dirs:      src
-  default-language:    Haskell2010
+  exposed-modules: FRP.Rhine.Bayes
+  other-modules: Data.Automaton.Bayes
+  build-depends:
+    automaton,
+    base >=4.14 && <4.20,
+    log-domain >=0.12,
+    mmorph ^>=1.2,
+    monad-bayes ^>=1.3,
+    rhine ==1.3,
+    transformers >=0.5
+
+  hs-source-dirs: src
+  default-language: Haskell2010
   default-extensions:
     Arrows
     DataKinds
@@ -51,24 +56,27 @@
     TypeFamilies
     TypeOperators
 
-  ghc-options:         -W
+  ghc-options: -W
+
   if flag(dev)
     ghc-options: -Werror
 
 executable rhine-bayes-gloss
-  main-is:             Main.hs
-  hs-source-dirs:      app
-  build-depends:       base         >= 4.11 && < 4.18
-                     , rhine
-                     , rhine-bayes
-                     , rhine-gloss  == 1.2
-                     , dunai
-                     , monad-bayes
-                     , transformers
-                     , log-domain
-                     , mmorph
-                     , time
-  default-language:    Haskell2010
+  main-is: Main.hs
+  hs-source-dirs: app
+  build-depends:
+    automaton,
+    base >=4.14 && <4.20,
+    log-domain,
+    mmorph,
+    monad-bayes,
+    rhine,
+    rhine-bayes,
+    rhine-gloss ==1.3,
+    time,
+    transformers
+
+  default-language: Haskell2010
   default-extensions:
     Arrows
     DataKinds
@@ -80,7 +88,12 @@
     TypeFamilies
     TypeOperators
 
-  ghc-options:         -W -threaded -rtsopts -with-rtsopts=-N
+  ghc-options:
+    -W
+    -threaded
+    -rtsopts
+    -with-rtsopts=-N
+
   if flag(dev)
     ghc-options: -Werror
 
diff --git a/src/Data/Automaton/Bayes.hs b/src/Data/Automaton/Bayes.hs
new file mode 100644
--- /dev/null
+++ b/src/Data/Automaton/Bayes.hs
@@ -0,0 +1,67 @@
+{-# LANGUAGE NamedFieldPuns #-}
+
+module Data.Automaton.Bayes where
+
+-- base
+import Control.Arrow
+
+-- transformers
+import Control.Monad.Trans.Reader (ReaderT (..))
+
+-- log-domain
+import Numeric.Log hiding (sum)
+
+-- monad-bayes
+import Control.Monad.Bayes.Population (PopulationT (..), fromWeightedList, runPopulationT)
+
+-- mmorph
+import Control.Monad.Morph (hoist)
+
+-- automaton
+import Data.Automaton (Automaton (..), handleAutomaton)
+import Data.Stream (StreamT (..))
+import Data.Stream.Result (Result (..))
+
+-- | Run the Sequential Monte Carlo algorithm continuously on an 'Automaton'
+runPopulationS ::
+  forall m a b.
+  (Monad m) =>
+  -- | Number of particles
+  Int ->
+  -- | Resampler
+  (forall x. PopulationT m x -> PopulationT m x) ->
+  Automaton (PopulationT m) a b ->
+  -- FIXME Why not Automaton m a (PopulationT b)
+  Automaton m a [(b, Log Double)]
+runPopulationS nParticles resampler =
+  handleAutomaton
+    ( runPopulationStream
+        (commuteReaderPopulation . hoist resampler . commuteReaderPopulationBack)
+        . hoist commuteReaderPopulation
+    )
+  where
+    commuteReaderPopulation :: forall m r a. (Monad m) => ReaderT r (PopulationT m) a -> PopulationT (ReaderT r m) a
+    commuteReaderPopulation = fromWeightedList . ReaderT . fmap runPopulationT . runReaderT
+
+    commuteReaderPopulationBack :: forall m r a. (Monad m) => PopulationT (ReaderT r m) a -> ReaderT r (PopulationT m) a
+    commuteReaderPopulationBack = ReaderT . fmap fromWeightedList . runReaderT . runPopulationT
+
+    runPopulationStream ::
+      forall m b.
+      (Monad m) =>
+      (forall x. PopulationT m x -> PopulationT m x) ->
+      StreamT (PopulationT m) b ->
+      StreamT m [(b, Log Double)]
+    runPopulationStream resampler StreamT {step, state} =
+      StreamT
+        { state = replicate nParticles (state, 1 / fromIntegral nParticles)
+        , step = \states -> do
+            resultsAndProbabilities <- runPopulationT $ normalize $ resampler $ do
+              state <- fromWeightedList $ pure states
+              step state
+            return $! Result (first resultState <$> resultsAndProbabilities) (first output <$> resultsAndProbabilities)
+        }
+
+-- FIXME see PR re-adding this to monad-bayes
+normalize :: (Monad m) => PopulationT m a -> PopulationT m a
+normalize = fromWeightedList . fmap (\particles -> second (/ (sum $ snd <$> particles)) <$> particles) . runPopulationT
diff --git a/src/Data/MonadicStreamFunction/Bayes.hs b/src/Data/MonadicStreamFunction/Bayes.hs
deleted file mode 100644
--- a/src/Data/MonadicStreamFunction/Bayes.hs
+++ /dev/null
@@ -1,53 +0,0 @@
-module Data.MonadicStreamFunction.Bayes where
-
--- base
-import Control.Arrow
-import Data.Functor (($>))
-import Data.Tuple (swap)
-
--- transformers
-
--- log-domain
-import Numeric.Log hiding (sum)
-
--- monad-bayes
-import Control.Monad.Bayes.Population
-
--- dunai
-import Data.MonadicStreamFunction
-import Data.MonadicStreamFunction.InternalCore (MSF (..))
-
--- | Run the Sequential Monte Carlo algorithm continuously on an 'MSF'
-runPopulationS ::
-  forall m a b.
-  (Monad m) =>
-  -- | Number of particles
-  Int ->
-  -- | Resampler
-  (forall x. PopulationT m x -> PopulationT m x) ->
-  MSF (PopulationT m) a b ->
-  -- FIXME Why not MSF m a (PopulationT b)
-  MSF m a [(b, Log Double)]
-runPopulationS nParticles resampler = runPopulationsS resampler . (spawn nParticles $>)
-
--- | Run the Sequential Monte Carlo algorithm continuously on a 'PopulationT' of 'MSF's
-runPopulationsS ::
-  (Monad m) =>
-  -- | Resampler
-  (forall x. PopulationT m x -> PopulationT m x) ->
-  PopulationT m (MSF (PopulationT m) a b) ->
-  MSF m a [(b, Log Double)]
-runPopulationsS resampler = go
-  where
-    go msfs = MSF $ \a -> do
-      -- TODO This is quite different than the dunai version now. Maybe it's right nevertheless.
-      -- FIXME This normalizes, which introduces bias, whatever that means
-      bAndMSFs <- runPopulationT $ normalize $ resampler $ flip unMSF a =<< msfs
-      return $
-        second (go . fromWeightedList . return) $
-          unzip $
-            (swap . fmap fst &&& swap . fmap snd) . swap <$> bAndMSFs
-
--- FIXME see PR re-adding this to monad-bayes
-normalize :: (Monad m) => PopulationT m a -> PopulationT m a
-normalize = fromWeightedList . fmap (\particles -> second (/ (sum $ snd <$> particles)) <$> particles) . runPopulationT
diff --git a/src/FRP/Rhine/Bayes.hs b/src/FRP/Rhine/Bayes.hs
--- a/src/FRP/Rhine/Bayes.hs
+++ b/src/FRP/Rhine/Bayes.hs
@@ -10,11 +10,11 @@
 import Control.Monad.Bayes.Class
 import Control.Monad.Bayes.Population
 
--- dunai
-import qualified Control.Monad.Trans.MSF.Reader as DunaiReader
+-- automaton
+import qualified Data.Automaton.Trans.Reader as AutomatonReader
 
--- dunai-bayes
-import qualified Data.MonadicStreamFunction.Bayes as DunaiBayes
+-- rhine-bayes
+import qualified Data.Automaton.Bayes as AutomatonBayes
 
 -- rhine
 import FRP.Rhine
@@ -24,18 +24,18 @@
 -- | Run the Sequential Monte Carlo algorithm continuously on a 'ClSF'.
 runPopulationCl ::
   forall m cl a b.
-  (Monad m) =>
+  (Monad m, MonadDistribution m) =>
   -- | Number of particles
   Int ->
   -- | Resampler (see 'Control.Monad.Bayes.PopulationT' for some standard choices)
-  (forall x. PopulationT m x -> PopulationT m x) ->
+  (forall x m. (MonadDistribution m) => PopulationT m x -> PopulationT m x) ->
   -- | A signal function modelling the stochastic process on which to perform inference.
   --   @a@ represents observations upon which the model should condition, using e.g. 'score'.
   --   It can also additionally contain hyperparameters.
   --   @b@ is the type of estimated current state.
   ClSF (PopulationT m) cl a b ->
   ClSF m cl a [(b, Log Double)]
-runPopulationCl nParticles resampler = DunaiReader.readerS . DunaiBayes.runPopulationS nParticles resampler . DunaiReader.runReaderS
+runPopulationCl nParticles resampler = AutomatonReader.readerS . AutomatonBayes.runPopulationS nParticles resampler . AutomatonReader.runReaderS
 
 -- * Short standard library of stochastic processes