0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

Typstで相関係数を算出

Posted at

Typstは配列のメソッドが結構充実しているので、ちょっとしたデータならRとか使わなくても基本的な記述統計とかはできるよね、ということで、ちょっとばかり遊んでみました。相関係数の練習問題(架空)とその解説の自動生成スクリプトです。

要素数5の配列(dataAdataB)の値を設定すると、問題用の表と、計算過程の説明つき解答を作成してくれます。解答の表示/非表示は、answer変数で切り替えます。

各種統計量の計算用スクリプト

平均値の算出

平均値は要素の総計を要素数で割れば良いので、配列のsumメソッドとlenメソッドの組み合わせで簡単に算出できます。

平均値の算出関数
#let calc_mean(array) = {array.sum() / array.len()}

標準偏差は分散の正の平方根で、(このデータ全体の)分散は平均値からの偏差の2乗値を平均したものなので、まずはmapメソッドで各要素の平均値からの偏差を求め、

偏差配列の生成関数
#let calc_deviation(array) = {array.map(x => x - calc_mean(array))}

それらを2乗した値の合計値を要素数で割って、その正の平方根を求めます。2乗してから合計する必要があるので、mapを使って要素の2乗値を求め、さらにsumで合計を求めてから、lenで割って平方根を求めています。

標準偏差の算出関数
#let calc_sd(array) = {
  calc.sqrt(
    calc_deviation(array).map(x => x * x ).sum()/array.len()
  )
}

なお、同じことはfoldでもできます。foldなら2乗の合計を一気に計算できますが、スクリプトの見た目はmapの方が若干短く、わかりやすいように思います。

標準偏差の算出関数 foldバージョン
#let calc_sd(array) = {
  calc.sqrt(
    calc_deviation(array).fold( 0, (sum, x) => sum + x*x )/array.len()
  )
}

なお、単に標準偏差を求めたいだけなら、この2つの処理を一度にまとめて行うこともできますが、共分散を求める際に両変数の偏差の配列が必要になるので、ここでは別々の関数にしてあります。

それから、相関係数(ピアソンの積率相関係数)はペアデータの共分散を標準化した値なので、共分散を求めるための関数も必要です。共分散は、ペア要素の(平均値からの)偏差の積の平均なので、まずは両変数の偏差の配列をzipでペアにしてから、mapを使ってペア要素の積を求め、sumで合計してから、lenで割って平均します。

共分散の算出関数
#let calc_cov(array1, array2) = {
  let paired = calc_deviation(array1).zip(calc_deviation(array2))
  paired.map(x => (x.at(0) * x.at(1))).sum() / paired.len()
}

標準偏差の場合と同様に、foldでも同じことができます。foldを使うとこうなります。

共分散の算出関数 foldバージョン
#let calc_cov(array1, array2) = {
  let paired = calc_deviation(array1).zip(calc_deviation(array2))
  paired.fold(0, (sum, x) => (sum + (x.at(0) * x.at(1)))) / paired.len()
}

そして、この共分散を両変数の標準偏差の積で割って標準化した値がピアソンの積率相関係数です。標準偏差を求める関数はすでに作成してありますし、たいした計算ではないのでわざわざ関数化しなくてもいいかもしれませんが、一応関数化しておきました。

相関係数の算出関数
#let calc_corr(array1, array2) = {
  calc_cov(array1, array2) / calc_sd(array1) / calc_sd(array2)
}

これで数値計算用関数の定義はおしまいです。

問題文の作成

あとは、typstの本文で問題文と解説文を作成するだけです。ここでは、問題文の直後にデータ用の配列(dataA、dataB)を作成し、それらの配列の要素と、先ほど作成した関数で求めた平均値、標準偏差をデータとして流し込んで表を作成しています。

typstの問題文
= 問題1

次のデータについて、ピアソンの積率相関係数 $r$ を求めよ。(値は小数第3位まで記入)

// データ
#let dataA = ( 5, 4, 9, 3, 7)
#let dataB = ( 7, 3, 5, 6, 3)

#figure()[
  #ztable(
    columns: (4em,4em,4em,4em,4em,4em,5em,5em),
    align: auto,
    toprule(),
    table.header(
      [変数], ..range(5).map( x => str(x+1)),[平均],[SD]
      ),
    midrule(),

    [A],
    ..dataA.map( x => str(x)), 
    num(digits:2)[#calc_mean(dataA)],
    num(digits:2)[#calc_sd(dataA)],

    [B],
    ..dataB.map( x => str(x)), 
    num(digits:2)[#calc_mean(dataB)],
    num(digits:2)[#calc_sd(dataB)],
    bottomrule()
  )
]

これで問題文はおしまいです。配列の値を変更すれば、自動で表の値も変更されます。

解説・解答欄の作成

ここからは解答欄および解説の部分です。ここではanswerというフラグ変数を定義して、answerの値で解説の表示・非表示を切り替えられるようにしました。こうしておけば、問題文と解説文のファイルが1つですみますし、問題の数値を訂正した場合に解説の数値が問題文とずれるようなこともなくなります。

ここから先は普通のtypstの文書なのでとくに解説はつけませんが、数式内で変数を使う場合、変数名にアンダースコア(_)が含まれているとトラブルの元なので、その点だけは注意が必要です。

typst
#let answer = true  // 解答表示あり: true、なし:false

#if answer [
  == 解答

  データAとデータBのピアソンの相関係数 $r$ は、次の式で定義される。

  $ "相関係数" r = frac("AとBの共分散","Aの標準偏差" times "Bの標準偏差") $

  また、共分散は次の式で定義される。

  $ "AとBの共分散" = frac(("A"_i"の偏差" times "B"_i"の偏差")"の合計", "データペアの個数") $

  #v(2em)

  上記の表のとおり、データAの平均値は#num(calc_mean(dataA), digits: 2)、データBの平均値は#num(calc_mean(dataB), digits: 2)であるから、各データの偏差、およびペアになるデータの偏差の積は次のとおりである。


  #let divA = calc_deviation(dataA)
  #let divB = calc_deviation(dataB)
  #let divAB = divA.zip(divB).map(x => x.at(0) * x.at(1))

  #figure()[
  //  #show table: format-table(auto, auto, auto)

    #ztable(
      columns: (6em,4em,4em,4em,4em,4em),
      format: (auto, (digits: 2), (digits: 2), (digits: 2), (digits: 2), (digits: 2),),
      align: horizon,
      toprule(),
      table.header(
        [変数], ..range(5).map( x => num(x+1))
        ),
      midrule(),

      [A],
      ..divA.map( x => str(x)), 

      [B],
      ..divB.map( x => str(x)), 

      midrule(),

      [偏差の積],
      ..divAB.map( x => str(x)),
      bottomrule()
    )
  ]

  したがって、

  $ "AとBの共分散" = frac( #divAB.map(x => num(x, digits:2)).join(" + ") , 5) = #num(calc_cov(dataA, dataB), digits: 2) $

  である。

  また、データAの標準偏差は#num(calc_sd(dataA), digits: 2)、データBの標準偏差は#num(calc_sd(dataB), digits: 2)であるから、相関係数 $r$ は次のように求められる。

  $ "相関係数" r = frac(
    #num(calc_cov(dataA, dataB), digits: 2), 
    #num(calc_sd(dataA), digits: 2) times 
    #num(calc_sd(dataB), digits: 2)) 

    = #num(calc_corr(dataA, dataB), digits: 4)dots.h
  $
]

#v(1fr)
#h(1fr)
#box(
  width: 8cm,
  stroke: (bottom: 1pt),
  inset: (bottom: 2pt),
  )[
    相関係数 $r$ :#h(1fr)#if answer [ #num(calc_corr(dataA, dataB), digits: 3)]#h(1fr)
  ]

まとめ

メソッドの練習がてらお遊び半分で作成してみましたが、Typstはなかなかに使い勝手がよいですね。これまでLaTeX + knitrで作成してきた文書の大部分がTypstだけで置き換えられそうです。

最後に、使用したスクリプト全文を掲載しておきます。

相関係数練習問題.typ
#import "@preview/zero:0.5.0": *
#import "@preview/booktabs:0.0.4": *
#show: booktabs-default-table-style

#let jfont = "Source Han Serif"
#set text(lang: "ja", font: jfont)

#show math.equation: set text(font: ("New Computer Modern Math", jfont))

///////////////////////////////////////////////////
// 計算用関数
///////////////////////////////////////////////////

// 平均値の算出関数
#let calc_mean(array) = {array.sum() / array.len()}

// 偏差配列の生成関数
#let calc_deviation(array) = {array.map(x => x - calc_mean(array))}

// 標準偏差の算出関数
#let calc_sd(array) = {
  calc.sqrt(
    calc_deviation(array).map(x => x * x ).sum()/array.len()
    // calc_deviation(array).fold( 0, (sum, x) => sum + x*x )/array.len() // こちらでも同じ
  )
}

// 共分散の算出関数
#let calc_cov(array1, array2) = {
  let paired = calc_deviation(array1).zip(calc_deviation(array2))
  paired.map(x => (x.at(0) * x.at(1))).sum() / paired.len()
  // paired.fold(0, (sum, x) => (sum + (x.at(0) * x.at(1)))) / paired.len() // こちらでも同じ
}



// 相関係数の算出関数
#let calc_corr(array1, array2) = {
  calc_cov(array1, array2) / calc_sd(array1) / calc_sd(array2)
}


///////////////////////////////////////////////////
// ここから本文
///////////////////////////////////////////////////
= 問題1

次のデータについて、ピアソンの積率相関係数 $r$ を求めよ。(値は小数第3位まで記入)

// データ
#let dataA = ( 5, 4, 9, 3, 7)
#let dataB = ( 7, 3, 5, 6, 3)

#figure()[
  #ztable(
    columns: (4em,4em,4em,4em,4em,4em,5em,5em),
    align: auto,
    toprule(),
    table.header(
      [変数], ..range(5).map( x => str(x+1)),[平均],[SD]
      ),
    midrule(),

    [A],
    ..dataA.map( x => str(x)), 
    num(digits:2)[#calc_mean(dataA)],
    num(digits:2)[#calc_sd(dataA)],

    [B],
    ..dataB.map( x => str(x)), 
    num(digits:2)[#calc_mean(dataB)],
    num(digits:2)[#calc_sd(dataB)],
    bottomrule()
  )
]

#v(2em)

#let answer = true  // 解答表示あり: true、なし:false

#if answer [
  == 解答

  データAとデータBのピアソンの相関係数 $r$ は、次の式で定義される。

  $ "相関係数" r = frac("AとBの共分散","Aの標準偏差" times "Bの標準偏差") $

  また、共分散は次の式で定義される。

  $ "AとBの共分散" = frac(("A"_i"の偏差" times "B"_i"の偏差")"の合計", "データペアの個数") $

  #v(2em)

  上記の表のとおり、データAの平均値は#num(calc_mean(dataA), digits: 2)、データBの平均値は#num(calc_mean(dataB), digits: 2)であるから、各データの偏差、およびペアになるデータの偏差の積は次のとおりである。


  #let divA = calc_deviation(dataA)
  #let divB = calc_deviation(dataB)
  #let divAB = divA.zip(divB).map(x => x.at(0) * x.at(1))

  #figure()[
  //  #show table: format-table(auto, auto, auto)

    #ztable(
      columns: (6em,4em,4em,4em,4em,4em),
      format: (auto, (digits: 2), (digits: 2), (digits: 2), (digits: 2), (digits: 2),),
      align: horizon,
      toprule(),
      table.header(
        [変数], ..range(5).map( x => num(x+1))
        ),
      midrule(),

      [A],
      ..divA.map( x => str(x)), 

      [B],
      ..divB.map( x => str(x)), 

      midrule(),

      [偏差の積],
      ..divAB.map( x => str(x)),
      bottomrule()
    )
  ]

  したがって、

  $ "AとBの共分散" = frac( #divAB.map(x => num(x, digits:2)).join(" + ") , 5) = #num(calc_cov(dataA, dataB), digits: 2) $

  である。

  また、データAの標準偏差は#num(calc_sd(dataA), digits: 2)、データBの標準偏差は#num(calc_sd(dataB), digits: 2)であるから、相関係数 $r$ は次のように求められる。

  $ "相関係数" r = frac(
    #num(calc_cov(dataA, dataB), digits: 2), 
    #num(calc_sd(dataA), digits: 2) times 
    #num(calc_sd(dataB), digits: 2)) 

    = #num(calc_corr(dataA, dataB), digits: 4)dots.h
  $
]

#v(1fr)
#h(1fr)
#box(
  width: 8cm,
  stroke: (bottom: 1pt),
  inset: (bottom: 2pt),
  )[
    相関係数 $r$ :#h(1fr)#if answer [ #num(calc_corr(dataA, dataB), digits: 3)]#h(1fr)
  ]

#v(2em)
0
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?