hmatrixのSTVectorとVectorとの部分更新速度比較

ノイズの中に信号を埋め込んで、その信号が検出できるかどうかのテストを行うことを考える。その場合、データに信号を埋め込むことを行うので、データ列の一部分に変更することになる。そこで、データ列の一部分の変更を行う関数を書くために、データ型の更新速度を調べてみた。今回は
hmatrixのデータ型STVectorとVectorとの比較をした。
方法としては、Vectorの方は、subVectorを用いて更新する部分としない部分に分けて、更新後にjoinするという方法

injectedV = join [subVector 0 n v
                       , add (subVector n nw v) w
                       , subVector (n+nw-1) (nv -nw -n) v]

をとり、STVectorの方は、Vector型のデータ列をunsafeThawVectorでSTVector型に変換して、in-placeの更新を行う方法

addInjsig n v w = runSTVector $ do
  v' <- unsafeThawVector v
  mapM_ (\i -> addInjsigCore v' w (n+i) i) [0 .. nw-1]
  return v'
  where
    nw = dim w

をとった。

テストコードは以下である。

--- testbench_STVector.hs
import Control.DeepSeq (deepseq)
import Control.Monad.ST (ST)
import Data.Packed.ST
import Numeric.LinearAlgebra
import System.Random
import Data.Time


main = do
  let v = fromList $ take 100000 $ randomRs (-10, 10) $mkStdGen 1 :: Vector Double
      w = fromList $ take 1000   $ randomRs (-10, 10) $mkStdGen 2 :: Vector Double
      n = 5000
      nw = dim w
      nv = dim v

  v `deepseq` return ()
  w `deepseq` return ()
  nw `deepseq` return ()
  nv `deepseq` return ()


  t1 <- getCurrentTime
  let injectedV = join [subVector 0 n v
                       , add (subVector n nw v) w
                       , subVector (n+nw-1) (nv -nw -n) v]
  print $ injectedV @> 5500
  t2 <- getCurrentTime
  print $ diffUTCTime t2 t1

  t3 <- getCurrentTime
  let injectedST = addInjsig n v w
  print $ injectedST @> 5500
  t4 <- getCurrentTime
  print $ diffUTCTime t4 t3



addInjsig n v w = runSTVector $ do
  v' <- unsafeThawVector v
  mapM_ (\i -> addInjsigCore v' w (n+i) i) [0 .. nw-1]
  return v'
  where
    nw = dim w

addInjsigCore :: STVector s Double -> Vector Double -> Int -> Int -> ST s ()
addInjsigCore v w i j = modifyVector v i (+w@>j)

結果は以下のようになった。

Vectorの方のコードは改善の余地がありそうだが、STVectorを用いたものはかなり速いという結果だった。

他の良いやり方のアイディアなどがあれば教えていただければありがたいです！