diff --git a/ChangeLog.md b/ChangeLog.md
--- a/ChangeLog.md
+++ b/ChangeLog.md
@@ -5,6 +5,11 @@
 ## Unreleased changes
 
 
+## Version 0.7.2.0
+
+-   `slynx`: Allow global normalization of mixture models.
+
+
 ## Version 0.7.1.0
 
 -   Be less strict with quoted identifiers/names in phylogenetic trees.
diff --git a/README.md b/README.md
--- a/README.md
+++ b/README.md
@@ -2,7 +2,7 @@
 
 # The ELynx Suite
 
-Version: 0.7.0.0.
+Version: 0.7.2.0.
 Reproducible evolution made easy.
 
 <p align="center"><img src="https://travis-ci.org/dschrempf/elynx.svg?branch=master"/></p>
@@ -73,13 +73,13 @@
     # OR: stack exec slynx -- --help
     # OR: slynx --help
 
-    Up to date
-    ELynx Suite version 0.7.0.0.
+    ELynx Suite version 0.7.2.0.
     Developed by Dominik Schrempf.
-    Compiled on July 10, 2022, at 08:48 am, UTC.
+    Compiled on May 3, 2023, at 13:18 pm, UTC.
     
     Usage: slynx [-v|--verbosity VALUE] [-o|--output-file-basename NAME]
                  [-f|--force] [--no-elynx-file] COMMAND
+    
       Analyze, and simulate multi sequence alignments.
     
     Available options:
diff --git a/elynx-markov.cabal b/elynx-markov.cabal
--- a/elynx-markov.cabal
+++ b/elynx-markov.cabal
@@ -1,6 +1,6 @@
 cabal-version:      3.0
 name:               elynx-markov
-version:            0.7.1.0
+version:            0.7.2.0
 synopsis:           Simulate molecular sequences along trees
 description:
   Examine, modify, and simulate molecular sequences in a reproducible way. Please see the README on GitHub at <https://github.com/dschrempf/elynx>.
diff --git a/src/ELynx/MarkovProcess/AminoAcid.hs b/src/ELynx/MarkovProcess/AminoAcid.hs
--- a/src/ELynx/MarkovProcess/AminoAcid.hs
+++ b/src/ELynx/MarkovProcess/AminoAcid.hs
@@ -368,12 +368,12 @@
 lgStatDist = pamlToAlphaVec lgStatDistPaml
 
 -- | LG substitution model.
-lg :: SubstitutionModel
-lg = substitutionModel Protein "LG" [] lgStatDist lgExch
+lg :: Normalize -> SubstitutionModel
+lg nz = substitutionModel Protein "LG" [] nz lgStatDist lgExch
 
 -- | LG substitution model with maybe a name and a custom stationary distribution.
-lgCustom :: Maybe String -> StationaryDistribution -> SubstitutionModel
-lgCustom mnm d = substitutionModel Protein nm [] d lgExch
+lgCustom :: Maybe String -> Normalize -> StationaryDistribution -> SubstitutionModel
+lgCustom mnm nz d = substitutionModel Protein nm [] nz d lgExch
   where
     nm = fromMaybe "LG-Custom" mnm
 
@@ -608,12 +608,12 @@
 wagStatDist = pamlToAlphaVec wagStatDistPaml
 
 -- | LG substitution model.
-wag :: SubstitutionModel
-wag = substitutionModel Protein "WAG" [] wagStatDist wagExch
+wag :: Normalize -> SubstitutionModel
+wag nz = substitutionModel Protein "WAG" [] nz wagStatDist wagExch
 
 -- | LG substitution model with maybe a name and a custom stationary distribution.
-wagCustom :: Maybe String -> StationaryDistribution -> SubstitutionModel
-wagCustom mnm d = substitutionModel Protein nm [] d wagExch
+wagCustom :: Maybe String -> Normalize -> StationaryDistribution -> SubstitutionModel
+wagCustom mnm nz d = substitutionModel Protein nm [] nz d wagExch
   where
     nm = fromMaybe "WAG-Custom" mnm
 
@@ -631,17 +631,17 @@
 uniformVec = V.replicate n (1 / fromIntegral n)
 
 -- | Poisson substitution model.
-poisson :: SubstitutionModel
-poisson = substitutionModel Protein "Poisson" [] uniformVec poissonExch
+poisson :: Normalize -> SubstitutionModel
+poisson nz = substitutionModel Protein "Poisson" [] nz uniformVec poissonExch
 
 -- | Poisson substitution model with maybe a name and a custom stationary distribution.
-poissonCustom :: Maybe String -> StationaryDistribution -> SubstitutionModel
-poissonCustom mnm d = substitutionModel Protein nm [] d poissonExch
+poissonCustom :: Maybe String -> Normalize -> StationaryDistribution -> SubstitutionModel
+poissonCustom mnm nz d = substitutionModel Protein nm [] nz d poissonExch
   where
     nm = fromMaybe "Poisson-Custom" mnm
 
 -- | General time reversible (GTR) substitution model for amino acids.
-gtr20 :: [Double] -> StationaryDistribution -> SubstitutionModel
-gtr20 es d = substitutionModel Protein "GTR" es d e
+gtr20 :: Normalize -> [Double] -> StationaryDistribution -> SubstitutionModel
+gtr20 nz es d = substitutionModel Protein "GTR" es nz d e
   where
     e = exchFromListUpper n es
diff --git a/src/ELynx/MarkovProcess/CXXModels.hs b/src/ELynx/MarkovProcess/CXXModels.hs
--- a/src/ELynx/MarkovProcess/CXXModels.hs
+++ b/src/ELynx/MarkovProcess/CXXModels.hs
@@ -117,6 +117,7 @@
 cxxSubstitutionModelFromStatDist nComps comp d =
   poissonCustom
     (Just name)
+    S.DoNotNormalize
     (normalizeSD d)
   where
     name = componentName nComps comp
@@ -139,6 +140,7 @@
 cxxFromStatDistsAndWeights ws ds =
   M.fromSubstitutionModels
     (cxxName n)
+    S.DoNormalize
     (V.fromList ws)
     sms
   where
diff --git a/src/ELynx/MarkovProcess/GammaRateHeterogeneity.hs b/src/ELynx/MarkovProcess/GammaRateHeterogeneity.hs
--- a/src/ELynx/MarkovProcess/GammaRateHeterogeneity.hs
+++ b/src/ELynx/MarkovProcess/GammaRateHeterogeneity.hs
@@ -58,8 +58,7 @@
     ++ "shape parameter "
     ++ show alpha
 
-splitSubstitutionModel ::
-  Int -> Double -> S.SubstitutionModel -> V.Vector S.SubstitutionModel
+splitSubstitutionModel :: Int -> Double -> S.SubstitutionModel -> V.Vector S.SubstitutionModel
 splitSubstitutionModel n alpha sm = renamedSMs
   where
     means = getMeans n alpha
@@ -67,9 +66,8 @@
     names = V.fromList $ map (("; gamma rate category " ++) . show) [1 :: Int ..]
     renamedSMs = V.zipWith S.appendName names scaledSMs
 
-expandSubstitutionModel ::
-  Int -> Double -> S.SubstitutionModel -> M.MixtureModel
-expandSubstitutionModel n alpha sm = M.fromSubstitutionModels name ws sms
+expandSubstitutionModel :: Int -> Double -> S.SubstitutionModel -> M.MixtureModel
+expandSubstitutionModel n alpha sm = M.fromSubstitutionModels name S.DoNotNormalize ws sms
   where
     name = S.name sm <> getName n alpha
     ws = V.replicate n 1.0
diff --git a/src/ELynx/MarkovProcess/MixtureModel.hs b/src/ELynx/MarkovProcess/MixtureModel.hs
--- a/src/ELynx/MarkovProcess/MixtureModel.hs
+++ b/src/ELynx/MarkovProcess/MixtureModel.hs
@@ -64,15 +64,33 @@
 getSubstitutionModels :: MixtureModel -> V.Vector S.SubstitutionModel
 getSubstitutionModels = V.map substModel . components
 
+normalizeGlobally :: V.Vector Weight -> V.Vector S.SubstitutionModel -> V.Vector S.SubstitutionModel
+normalizeGlobally ws ss = V.map (S.scale $ recip c) ss
+  where
+    cks = V.map S.totalRate ss
+    cNoWeights = V.sum $ V.zipWith (*) ws cks
+    c = cNoWeights / V.sum ws
+
 -- | Create a mixture model from a list of substitution models.
-fromSubstitutionModels :: S.Name -> V.Vector Weight -> V.Vector S.SubstitutionModel -> MixtureModel
-fromSubstitutionModels n ws sms
+--
+-- If 'S.Normalize' is 'S.DoNormalize', globally normalize the mixture model.
+-- Global normalization has no effect if all components are already normalized.
+fromSubstitutionModels ::
+  S.Name ->
+  S.Normalize ->
+  V.Vector Weight ->
+  V.Vector S.SubstitutionModel ->
+  MixtureModel
+fromSubstitutionModels n nz ws sms
   | null ws = error "fromSubstitutionModels: No weights given."
   | length ws /= length sms = error "fromSubstitutionModels: Number of weights and substitution models does not match."
   | not $ allEqual alphs = error "fromSubstitutionModels: alphabets of substitution models are not equal."
   | otherwise = MixtureModel n (V.head alphs) comps
   where
-    comps = V.zipWith Component ws sms
+    smsNormalized = case nz of
+      S.DoNormalize -> normalizeGlobally ws sms
+      S.DoNotNormalize -> sms
+    comps = V.zipWith Component ws smsNormalized
     alphs = V.map S.alphabet sms
     allEqual xs
       | V.null xs = True
@@ -80,7 +98,7 @@
 
 -- | Concatenate mixture models.
 concatenate :: S.Name -> V.Vector MixtureModel -> MixtureModel
-concatenate n mms = fromSubstitutionModels n ws sms
+concatenate n mms = fromSubstitutionModels n S.DoNotNormalize ws sms
   where
     comps = V.concatMap components mms
     ws = V.map weight comps
diff --git a/src/ELynx/MarkovProcess/Nucleotide.hs b/src/ELynx/MarkovProcess/Nucleotide.hs
--- a/src/ELynx/MarkovProcess/Nucleotide.hs
+++ b/src/ELynx/MarkovProcess/Nucleotide.hs
@@ -81,12 +81,12 @@
 uniformVec = V.replicate n (1 / fromIntegral n)
 
 -- | JC substitution model.
-jc :: SubstitutionModel
-jc = substitutionModel DNA "JC" [] uniformVec jcExch
+jc :: Normalize -> SubstitutionModel
+jc nz = substitutionModel DNA "JC" [] nz uniformVec jcExch
 
 -- | F81 substitution model.
-f81 :: StationaryDistribution -> SubstitutionModel
-f81 d = substitutionModel DNA "F81" [] d jcExch
+f81 :: Normalize -> StationaryDistribution -> SubstitutionModel
+f81 nz d = substitutionModel DNA "F81" [] nz d jcExch
 
 hkyExch :: Double -> ExchangeabilityMatrix
 hkyExch k =
@@ -94,11 +94,11 @@
     [0.0, 1.0, k, 1.0, 1.0, 0.0, 1.0, k, k, 1.0, 0.0, 1.0, 1.0, k, 1.0, 0.0]
 
 -- | HKY substitution model.
-hky :: Double -> StationaryDistribution -> SubstitutionModel
-hky k d = substitutionModel DNA "HKY" [k] d e where e = hkyExch k
+hky :: Normalize -> Double -> StationaryDistribution -> SubstitutionModel
+hky nz k d = substitutionModel DNA "HKY" [k] nz d e where e = hkyExch k
 
 -- | HKY substitution model.
-gtr4 :: [Double] -> StationaryDistribution -> SubstitutionModel
-gtr4 es d = substitutionModel DNA "GTR" es d e
+gtr4 :: Normalize -> [Double] -> StationaryDistribution -> SubstitutionModel
+gtr4 nz es d = substitutionModel DNA "GTR" es nz d e
   where
     e = exchFromListUpper n es
diff --git a/src/ELynx/MarkovProcess/SubstitutionModel.hs b/src/ELynx/MarkovProcess/SubstitutionModel.hs
--- a/src/ELynx/MarkovProcess/SubstitutionModel.hs
+++ b/src/ELynx/MarkovProcess/SubstitutionModel.hs
@@ -27,6 +27,7 @@
     totalRate,
 
     -- * Building substitution models
+    Normalize (..),
     substitutionModel,
 
     -- * Transformations
@@ -92,17 +93,26 @@
     s = V.sum v
 {-# INLINE normalizeSumVec #-}
 
+-- | Should the substitution model be normalized?
+data Normalize = DoNormalize | DoNotNormalize
+  deriving (Show, Read, Eq)
+
 -- | Create normalized 'SubstitutionModel'. See 'normalize'.
 substitutionModel ::
   Alphabet ->
   Name ->
   Params ->
+  Normalize ->
   R.StationaryDistribution ->
   R.ExchangeabilityMatrix ->
   SubstitutionModel
-substitutionModel c n ps d e =
+substitutionModel c n ps nz d e =
   if R.isValid d
-    then normalize $ SubstitutionModel c n ps (normalizeSumVec d) e
+    then
+      let sm = SubstitutionModel c n ps (normalizeSumVec d) e
+       in case nz of
+            DoNormalize -> normalize sm
+            DoNotNormalize -> sm
     else
       error $
         "substitionModel: Stationary distribution does not sum to 1.0: "
diff --git a/test/ELynx/MarkovProcess/AminoAcidSpec.hs b/test/ELynx/MarkovProcess/AminoAcidSpec.hs
--- a/test/ELynx/MarkovProcess/AminoAcidSpec.hs
+++ b/test/ELynx/MarkovProcess/AminoAcidSpec.hs
@@ -16,6 +16,7 @@
 import qualified Data.Vector.Storable as V
 import ELynx.MarkovProcess.AminoAcid
 import qualified ELynx.MarkovProcess.RateMatrix as R
+import ELynx.MarkovProcess.SubstitutionModel (Normalize (DoNormalize))
 import qualified ELynx.MarkovProcess.SubstitutionModel as S
 import ELynx.Tools.Equality
 import Numeric.LinearAlgebra
@@ -482,13 +483,13 @@
 statDistUniform = vector $ replicate 20 0.05
 
 statDistLG :: R.StationaryDistribution
-statDistLG = S.stationaryDistribution lg
+statDistLG = S.stationaryDistribution $ lg DoNormalize
 
 exchLG :: R.ExchangeabilityMatrix
-exchLG = S.exchangeabilityMatrix lg
+exchLG = S.exchangeabilityMatrix $ lg DoNormalize
 
 rmLG :: R.RateMatrix
-rmLG = S.rateMatrix lg
+rmLG = S.rateMatrix $ lg DoNormalize
 
 spec :: Spec
 spec = do
@@ -513,11 +514,12 @@
               S.rateMatrix $
                 lgCustom
                   Nothing
+                  DoNormalize
                   statDistUniform
       f `nearlyEqVec` statDistUniform `shouldBe` True
   describe "poisson" $
     it "stationary distribution is uniform 1/20" $
-      R.getStationaryDistribution (S.rateMatrix poisson)
+      R.getStationaryDistribution (S.rateMatrix $ poisson DoNormalize)
         `nearlyEqVec` statDistUniform
         `shouldBe` True
   describe "poissonCustom" $
@@ -527,5 +529,6 @@
               S.rateMatrix $
                 poissonCustom
                   Nothing
+                  DoNormalize
                   statDistLGPython
       f `nearlyEqVec` statDistLGPython `shouldBe` True
diff --git a/test/ELynx/MarkovProcess/NucleotideSpec.hs b/test/ELynx/MarkovProcess/NucleotideSpec.hs
--- a/test/ELynx/MarkovProcess/NucleotideSpec.hs
+++ b/test/ELynx/MarkovProcess/NucleotideSpec.hs
@@ -24,7 +24,7 @@
 stationaryDist = fromList [0.2, 0.3, 0.3, 0.2]
 
 hkyModel :: SubstitutionModel
-hkyModel = hky 6.0 stationaryDist
+hkyModel = hky DoNormalize 6.0 stationaryDist
 
 spec :: Spec
 spec =
diff --git a/test/ELynx/Simulate/MarkovProcessAlongTreeSpec.hs b/test/ELynx/Simulate/MarkovProcessAlongTreeSpec.hs
--- a/test/ELynx/Simulate/MarkovProcessAlongTreeSpec.hs
+++ b/test/ELynx/Simulate/MarkovProcessAlongTreeSpec.hs
@@ -15,6 +15,7 @@
 
 import Data.Tree
 import ELynx.MarkovProcess.Nucleotide
+import ELynx.MarkovProcess.SubstitutionModel (Normalize (DoNormalize))
 import qualified ELynx.MarkovProcess.SubstitutionModel as S
 import ELynx.Simulate.MarkovProcess
 import ELynx.Simulate.MarkovProcessAlongTree
@@ -50,5 +51,5 @@
     tr <- simulate 10 d e testTree gen
     (length . rootLabel $ tr) `shouldBe` 10
   where
-    d = S.stationaryDistribution jc
-    e = S.exchangeabilityMatrix jc
+    d = S.stationaryDistribution $ jc DoNormalize
+    e = S.exchangeabilityMatrix $ jc DoNormalize