出力(プログラムで処理した結果)を「データフレーム」に入れる

データの種類

変数:箱に一つだけ値が入る

配列:複数の箱が並んだ変数、Rでは1次元の配列を「ベクトル」と呼ぶ

データフレーム:配列をまとめたデータセット

配列とデータフレームの図

Type, Token, TTR, GI, AWL, ASLに文の数(NoS)を加えて7つの言語指標をdata.frameに

  • Number of Sentences は、data.tmpの要素数

流れ:各指標のベクトルを作っておいて、できたものをデータフレームにまとめる: data.frame()

  1. ベクトルの入れ物を用意(初期化)
  2. それぞれのベクトルに数値を入れていく
  3. ベクトルをまとめてデータフレームにする
myIndexesDf <- function(){
  
  # 各指標のベクトルを作っておいて、できたものをデータフレームにまとめる
  
  # ベクトルの初期化
  typeV  <- c()
  tokenV <- c()
  ttrV   <- c()
  giV    <- c()
  awlV   <- c()
  aslV   <- c()
  nosV   <- c()
  
  files <- list.files()
  
  for (i in files) {
    lines.tmp <- scan(i, what="char", sep="\n")
    body.tmp <- grep("^\\*\\w+:\t", lines.tmp, value=T)
    data.tmp <- gsub("^\\*\\w+:\t", "", body.tmp)
    data.tmp <- data.tmp[data.tmp != ""]
    tmp4 <- gsub("[[:punct:]]", " ", data.tmp) 
    tmp5 <- gsub("  +", " ", tmp4)
    tmp6 <- tolower(tmp5)
    tmp7 <- strsplit(tmp6, " ")
    tmp8 <- unlist(tmp7)
    token <- sort(tmp8)
    type <- unique(token)
    
    ttr <- length(type)/length(token)
    gi  <- length(type)/sqrt(length(token))
    tmp   <- paste(token, collapse = "")
    awl <- nchar(tmp)/length(token)
    asl <- length(token)/length(data.tmp)
    nos <- length(body.tmp)                                # Number of Sentences
    
#   cat(length(type), length(token), ttr, gi, awl, asl,"\n")
    
    # 各数値をベクトルに入れる
    typeV <- c(typeV, length(type))
    tokenV <- c(tokenV, length(token))
    ttrV   <- c(ttrV, ttr)
    giV    <- c(giV, gi)
    awlV   <- c(awlV, awl)
    aslV   <- c(aslV, asl)
    nosV   <- c(nosV, nos)
  

  }
  
 # 最後に全部をデータフレームにまとめる
 data.frame(typeV, tokenV, ttrV, giV, awlV, aslV, nosV)

}

実行

setwd("NICER1_3/NICER_NS")

myIndexesDf()

出来上がったデータフレームをオブジェクトとして「保存」

setwd("NICER1_3/NICER_NS")

indexes.df <- myIndexesDf()

indexes.df
head(indexes.df)

「見出し」を付け替える:names() <- c()

names(indexes.df) <- c("Type","Token","TTR","GI","AWL","ASL", "NoS")

head(indexes.df)

データフレームの基本操作

データフレームの構造の確認:str()

str(indexes.df)
## 'data.frame':    71 obs. of  7 variables:
##  $ Type : int  359 343 356 340 396 337 264 334 370 300 ...
##  $ Token: int  742 648 849 832 908 845 605 852 1092 763 ...
##  $ TTR  : num  0.484 0.529 0.419 0.409 0.436 ...
##  $ GI   : num  13.2 13.5 12.2 11.8 13.1 ...
##  $ AWL  : num  4.56 5.1 5.47 5.23 4.7 ...
##  $ ASL  : num  19 24.9 38.6 27.7 23.3 ...
##  $ NoS  : int  39 26 22 30 39 31 27 43 57 22 ...

データフレームの一部(特定の列(カラム))だけ取り出す: データフレーム$見出し

indexes.df$Token
##  [1]  742  648  849  832  908  845  605  852 1092  763  887  604  531  916 1154
## [16]  904  672 1158 2340 1133  965 1078  944 1094 1060 1496  977  974  805 1019
## [31]  943 1073 1346  893 1545  800 1080 1059 1091 1482 1024 1050  736  512  676
## [46] 1187 1197  559 1503  547 1163  456 1822  701 1037  891 1088  963 1187  714
## [61]  786  991  962 1181  966 1229  557  947 1133  965 1156

ファイルに保存:write.table(データフレーム, “ファイル名”)

write.table(indexes.df, "indexes.df.txt", row.names=F)
  • 区切り文字はスペース
    • カンマに変える場合は、オプション sep=“,”
  • 左に通し番号が入る
    • 入れないオプションは、row.names=F

基本統計量を見る:summary()

summary(indexes.df)
##       Type           Token             TTR               GI        
##  Min.   :222.0   Min.   : 456.0   Min.   :0.3116   Min.   : 9.928  
##  1st Qu.:336.5   1st Qu.: 802.5   1st Qu.:0.3608   1st Qu.:11.507  
##  Median :383.0   Median : 965.0   Median :0.4037   Median :12.288  
##  Mean   :387.3   Mean   : 986.5   Mean   :0.4044   Mean   :12.376  
##  3rd Qu.:428.0   3rd Qu.:1113.5   3rd Qu.:0.4362   3rd Qu.:13.269  
##  Max.   :761.0   Max.   :2340.0   Max.   :0.5293   Max.   :15.732  
##       AWL             ASL             NoS        
##  Min.   :4.013   Min.   :11.75   Min.   : 22.00  
##  1st Qu.:4.572   1st Qu.:19.10   1st Qu.: 34.50  
##  Median :4.774   Median :21.93   Median : 43.00  
##  Mean   :4.751   Mean   :22.59   Mean   : 45.66  
##  3rd Qu.:4.900   3rd Qu.:25.00   3rd Qu.: 50.50  
##  Max.   :5.571   Max.   :38.59   Max.   :155.00

一部分だけの基本統計量を見る:Tokenを例に

summary(indexes.df$Token)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   456.0   802.5   965.0   986.5  1113.5  2340.0

基本的なグラフ

boxplot()

• 真ん中の太い線:中央値 Median(50%)、平均ではないので注意
• 箱の範囲:25%から75%の範囲(これを四分範囲と呼ぶ)
    (つまり、半分のデータはこの範囲に入る)
• ひげの範囲:箱の端から箱の長さの1.5倍以内にある「実際の」数値の最大値と最小値
• 外れている○印:外れ値(ひげの範囲外のもの)
boxplot(indexes.df)

例:個別のグラフ:TypeとTokenの箱ひげ図

boxplot(indexes.df$Type, indexes.df$Token)

散布図: plot()

  • すべての組み合わせの2次元の散布図
plot(indexes.df)

例:個別のグラフ:TypeとTokenの散布図

plot(indexes.df$Type, indexes.df$Token)

棒グラフ: barplot()

  • 基本は一つのベクトルだけ
barplot(indexes.df$Token)

例:並べ替えて棒グラフにする

barplot(sort(indexes.df$Token))

ヒストグラム: hist()

  • 基本は一つのベクトルだけ
hist(indexes.df$Token)

♪ 演習

演習1:TokenとTTRの関係を散布図に

演習2:TokenとGIの関係を散布図に


解答例1

plot(indexes.df$Token, indexes.df$TTR)

回帰式を描いてみる: abline(lm(y~x))

  • abline() 直線を書く(a地点からb地点へ?)
  • lm(y~x) Linear Model(線形モデル) 目的変数 ~ 説明変数 の関係
    • どっちが「目的」(Y軸)でどっちが「説明」(X軸)か考えて
plot(indexes.df$Token, indexes.df$TTR)
abline(lm(indexes.df$TTR~indexes.df$Token))

plot(indexes.df$Token, indexes.df$GI)
abline(lm(indexes.df$GI~indexes.df$Token))

二つを合わせて描いてみる

  • par(new=T) で二つ目を上に重ね描きする
  • col=色オプションで、色を変える
plot(indexes.df$Token, indexes.df$TTR)
abline(lm(indexes.df$TTR~indexes.df$Token))

par(new=T)

plot(indexes.df$Token, indexes.df$GI, col="red")
abline(lm(indexes.df$GI~indexes.df$Token), col="red")

  • X軸は同じ、Y軸が違う点に注意

もう一つの語彙多様性指標:MATTR (Moving-Average Type-Token Ratio)

TTR, GIの問題点と解決案

  • テキストの長さ(token)の影響を受ける
  • テキストの長さを一定にして、比較すればよい
  • テキストの長さは必ずしも同じではない
  • 特定の長さで切って、平均を出せばよい
  • 特定の長さを統一にすれば、テキスト間で比較できる

Covington and McFall (2008, 2010)

Michael A. Covington & Joe D. McFall (2010) Cutting the Gordian Knot: The Moving-Average Type–Token Ratio (MATTR), Journal of Quantitative Linguistics, 17:2, 94-100, DOI: 10.1080/09296171003643098

Covington and McFall (2008) MATTR

MATTR

# 阿部大輔 2014, 杉浦正利 2020 NICER用に話者記号の削除・修正
mattr <- function () {
    wttr <- 0
    ttrsum <- 0

    lines.tmp <- scan(choose.files(), what="char", sep="\n")
    # body部分の抽出と話者記号と句読点の削除
  body.tmp <- grep("^\\*\\w+:\t", lines.tmp, value=T)
  data.tmp <- gsub("^\\*\\w+:\t", "", body.tmp)
  data.tmp <- data.tmp[data.tmp != ""]
  tmp4 <- gsub("[[:punct:]]", " ", data.tmp) 
  tmp5 <- gsub("  +", " ", tmp4)
    
    token <- tolower(unlist(strsplit(tmp5, "\\W+")))
    
    numwords <- length(token)                              # 元々の語数を記録しておく
    
    token <- c(token,token)                                        # yを一周するとまた先頭の語が出てくるようにyとyをつなげる
    
    for(i in 1:numwords){                             #↑で記録した総語数の回数ループする
        mado <- token[i:(99+i)]                           # 100語の幅で順に右にずらしていく
        wttr <- length(unique(sort(mado)))/100           #100語ずつ=tokenが100
        ttrsum <- ttrsum + wttr                          # 出てきたTTRを足していく
    }
    ttrsum/numwords                                   # 総TTR数を総語数の回数で割って平均を出す
}
#setwd("NICER1_3/NICER_NS")

mattr()

♪ 演習

演習3:MATTRを修正して、フォルダー内のファイルすべてについてMATTRが出るようにする。

  • mattr2() という名前でスクリプト作成
mattr2 <- function () {
  
  mattrV <- c()                                   # ベクトルデータの準備
  
    wttr <- 0
    ttrsum <- 0

    files <- list.files()                            # ファイルのリスト
  
  for (i in files) {                               # すべてのファイルに対して
    
    lines.tmp <- scan(i, what="char", sep="\n")
  body.tmp <- grep("^\\*\\w+:\t", lines.tmp, value=T)
  data.tmp <- gsub("^\\*\\w+:\t", "", body.tmp)
  data.tmp <- data.tmp[data.tmp != ""]
  tmp4 <- gsub("[[:punct:]]", " ", data.tmp) 
  tmp5 <- gsub("  +", " ", tmp4)
    
    token <- tolower(unlist(strsplit(tmp5, "\\W+")))
    
    numwords <- length(token)                              # 元々の語数を記録しておく
    
    token <- c(token,token)                                        # yを一周するとまた先頭の語が出てくるようにyとyをつなげる
    
    for(i in 1:numwords){                             #↑で記録した総語数の回数ループする
        mado <- token[i:(99+i)]                           # 100語の幅で順に右にずらしていく
        wttr <- length(unique(sort(mado)))/100           #100語ずつ=tokenが100
        ttrsum <- ttrsum + wttr                          # 出てきたTTRを足していく
    }
    #ttrsum/numwords                                   # 総TTR数を総語数の回数で割って平均を出す
    mattr <-    ttrsum/numwords                         # mattrに代入

    mattrV   <- c(mattrV, mattr)                      # ベクトルデータとして保存
  }                                                 # forに対応して閉じる
mattrV                                              # ベクトルデータを表示・出力
}
#setwd("NICER1_3/NICER_NS")    # 母語話者データ
#setwd("NICER1_3/NICER_NNS")   # 学習者データ
mattr2()

データフレームへのデータの追加:data.frame(もとのdf, 追加するベクトルデータ)

例:mattr2() の結果のベクトルデータを、データフレームidexes.dfに追加する

setwd("NICER1_3/NICER_NS")

mattr <- mattr2()

indexes.df2 <- data.frame(indexes.df, mattr)          # 

head(indexes.df2)

演習4:MATTRはテキストの長さの影響を受けないことを散布図を作って確認する

plot(indexes.df2$Token, indexes.df2$mattr)
abline(lm(indexes.df2$mattr~indexes.df2$Token))

plot(indexes.df$Token, indexes.df$TTR)
abline(lm(indexes.df$TTR~indexes.df$Token))

par(new=T)

plot(indexes.df$Token, indexes.df$GI, col="red")
abline(lm(indexes.df$GI~indexes.df$Token), col="red")

par(new=T)

plot(indexes.df2$Token, indexes.df2$mattr, col="blue")
abline(lm(indexes.df2$mattr~indexes.df2$Token), col="blue")

  • TTRやGIよりもMATTRの線の傾きが緩やかなので、改善されているといえるが、まだ傾いている。
  • データの性質による可能性がある。
    • 別々のデータで多様性指標を出しているので、そもそも多様性が同じではないかも。
    • ということは、多様性が等しいはずのデータをもとに、テキスト長を変えて、指標の変化を見る必要がある。

ヘッダー情報を取り込む:学習者要因やエッセイのスコアとの関係の分析

ヘッダー情報の取り込み

データフォーマットを考える

  • ヘッダー部分は、行頭が @ で始まっている

学習者データのヘッダー部分の例

@Begin
@Participants:  JPN501
@PID:   PIDJP501
@Age:   21
@Sex:   F
@YearInSchool:  U2
@Major: agriculture
@StudyHistory:  8
@OtherLanguage: Chinese=1.0;none=
@Qualification: TOEIC=590(2013);none=;none=
@Abroad:    none=;none=
@Reading:   3
@Writing:   2
@Listening: 2
@Speaking:  1
@JapaneseEssay: 4
@EnglishEssayEx:    3
@EnglishEssay:  2
@Difficulty:    
@EssayTraining: 3
@SelfEval:  2
@TopicEase: 4
@Topic: sports
@Criterion: 4
@Proctor:   1
@Comments:  
@Date:  2013-12-17
@Version:   1.3

母語話者データのヘッダー部分の例

@Begin
@Participants:  NS501
@PID:   PIDNS501
@Age:   27
@Sex:   M
@L1:    AmE
@FatherL1:  none
@MotherL1:  none
@AcademicBackground:    M1
@OtherLanguage: Japanese=0.7=good;none==
@Topic: education
@EnglishEssay:  5
@SelfEval:  1
@TopicEase: 3
@EssayTraining: 5
@Proctor:   1
@Criterion: 5
@Comments:  
@Date:  2013-12-02
@Version:   1.3

エッセイ評価(Criterionのスコア)を取り込んでみる

考え方

  1. 検索(lines.tmpに一文一行でデータが入っているとして)
criterion.tmp <- grep("@Criterion", lines.tmp, value = T)
  1. 置換・削除:その行の値のところだけ残して前を消す
score <- gsub("@Criterion:\t", "", criterion.tmp)

♪ 演習

演習5:フォルダー内のファイルすべてについてスコアが出るようにする。

  • myScore() という名前でスクリプト作成

(^^♪ 解答例

setwd("NICER1_3/NICER_NS")

myScore <- function(){
  
  scoreV <- c()                                             # ベクトルscoreVの初期化
  
  files <- list.files()
  for (i in files) {

  lines.tmp <- scan(i, what="char", sep="\n")                  

  criterion.tmp <- grep("@Criterion", lines.tmp, value = T) # @Criterionの行を検索
 
  score <- gsub("@Criterion:\t", "", criterion.tmp)         # スコアのみに

   scoreV   <- c(scoreV , score)                             # ベクトルデータとして保存
  }
 scoreV
}

実行

setwd("NICER1_3/NICER_NS")

myScore()

演習6:score() の結果を、ベクトルに保存し、それをデータフレームに追加する

setwd("NICER1_3/NICER_NS")

score <- myScore()


indexes.df3 <- data.frame(indexes.df2, score)         

head(indexes.df3)

演習7:参加者ID(JPN501とかNS501)も、データフレームに追加する

setwd("NICER1_3/NICER_NS")

myID <- function(){
  
  idV <- c()                                             # ベクトルidVの初期化
  
  files <- list.files()
  for (i in files) {

  lines.tmp <- scan(i, what="char", sep="\n")                  

  id.tmp <- grep("@Participants:", lines.tmp, value = T) # @Participantsの行を検索
 
  id <- gsub("@Participants:\t", "", id.tmp)         # IDのみに

   idV   <- c(idV , id)                             # ベクトルデータとして保存
  }
 idV
}
setwd("NICER1_3/NICER_NS")

id <- myID()


indexes.df4 <- data.frame(indexes.df3, id)          

head(indexes.df4)

演習8:エッセイのスコアと各種言語指標との関係を散布図を作って確認する

plot(indexes.df4)

演習9:ここまでのまとめ。IDと「トピック」とCriterionスコアと8つの言語指標を一度にまとめて出してくれるスクリプトを書いてみましょう。

  • それができれば、学習者データのデータフレームもすぐに出ます。
myNICERindex <- function(){
  
  # 各指標のベクトルを作っておいて、できたものをデータフレームにまとめる
  
  # ベクトルの初期化
  typeV  <- c()
  tokenV <- c()
  ttrV   <- c()
  giV    <- c()
  awlV   <- c()
  aslV   <- c()
  nosV   <- c()
  
  mattrV <- c()                                   # MATTR追加 以下同じ
  scoreV <- c()
  idV    <- c()
  topicV <- c()
  
  files <- list.files()
  
  for (i in files) {
    lines.tmp <- scan(i, what="char", sep="\n")
    body.tmp <- grep("^\\*\\w+:\t", lines.tmp, value=T)
    data.tmp <- gsub("^\\*\\w+:\t", "", body.tmp)
    data.tmp <- data.tmp[data.tmp != ""]
    tmp4 <- gsub("[[:punct:]]", " ", data.tmp) 
    tmp5 <- gsub("  +", " ", tmp4)
    tmp6 <- tolower(tmp5)
    tmp7 <- strsplit(tmp6, " ")
    tmp8 <- unlist(tmp7)
    token <- sort(tmp8)
    type <- unique(token)
    
    ttr <- length(type)/length(token)
    gi  <- length(type)/sqrt(length(token))
    tmp <- paste(token, collapse = "")
    awl <- nchar(tmp)/length(token)
    asl <- length(token)/length(data.tmp)
    nos <- length(body.tmp)                                # Number of Sentences
    
    #MATTR 追加
    numwords <- length(token)
    token <- c(token,token)
    ttrsum <- 0                                            # ttrsum の初期化
    mattr <- 0                                             # mattr の初期化
    for(i in 1:numwords){
      mado <- token[i:(99+i)]
      wttr <- length(unique(sort(mado)))/100 
      ttrsum <- ttrsum + wttr
    }
    mattr <-    ttrsum/numwords                         # mattrに代入
    mattrV   <- c(mattrV, mattr)                      # ベクトルデータとして保存
    
    #スコア 追加
    criterion.tmp <- grep("@Criterion", lines.tmp, value = T) # @Criterionの行を検索
    score <- gsub("@Criterion:\t", "", criterion.tmp)         # スコアのみに
    scoreV   <- c(scoreV , score)                             # ベクトルデータとして保存
    
    #ID 追加
    id.tmp <- grep("@Participants:", lines.tmp, value = T) # @Participantsの行を検索
    id <- gsub("@Participants:\t", "", id.tmp)         # IDのみに
    idV   <- c(idV , id)                             # ベクトルデータとして保存
 
    #Topic 追加
    topic.tmp <- grep("@Topic:", lines.tmp, value = T) # @Participantsの行を検索
    topic <- gsub("@Topic:\t", "", topic.tmp)         # IDのみに
    topicV   <- c(topicV , topic)                             # ベクトルデータとして保存
  
    
    # 各数値をベクトルに入れる
    typeV <- c(typeV, length(type))
    tokenV <- c(tokenV, length(token))
    ttrV   <- c(ttrV, ttr)
    giV    <- c(giV, gi)
    awlV   <- c(awlV, awl)
    aslV   <- c(aslV, asl)
    nosV   <- c(nosV, nos)
  

  }
  
 # 最後に全部をデータフレームにまとめる
 data.frame(idV, topicV, scoreV, typeV, tokenV, ttrV, giV, mattrV, awlV, aslV, nosV)

}

実行

母語話者データ

setwd("NICER1_3/NICER_NS")

ns_indexes.df <- myNICERindex()
head(ns_indexes.df )

「見出し」を付け替える:names() <- c()

names(ns_indexes.df ) <- c("ID", "Topic", "Score", "Type", "Token", "TTR", "GI", "MATTR", "AWL","ASL", "NoS")

head(ns_indexes.df )

学習者データ

setwd("NICER1_3/NICER_NNS")                    # ディレクトリーを変える NNS

jp_indexes.df <- myNICERindex()
names(jp_indexes.df ) <- c("ID", "Topic", "Score", "Type", "Token", "TTR", "GI", "MATTR", "AWL","ASL", "NoS")
head(jp_indexes.df )

学習者データと母語話者データの統合

学習者データのデータフレームと母語話者データのデータフレームの作成

  • jp_indexes.df と ns_indexes.df

nsとjpというカテゴリーをそれぞれのデータフレームに追加する

  • 統合した後、カテゴリーでデータを分けることができるように

  • 新しいカラムを追加:使ってない見出しを設定し、そこに入れたいものを代入

    • ns と jp は、言語の違いということで、L1とL2、Langという見出しにして、ns は 1 、jp は 2 を入れる。
jp_indexes.df$Lang <- 2

head(jp_indexes.df)
ns_indexes.df$Lang <- 1

head(ns_indexes.df)

二つのデータフレームを統合する:「縦に(row)」つなげる: rbind(データフレーム1, データフレーム2)

  • all_indexes.df
all_indexes.df <- rbind(jp_indexes.df, ns_indexes.df)

head(all_indexes.df)
tail(all_indexes.df)

データ構造の確認

str(all_indexes.df)
## 'data.frame':    1214 obs. of  12 variables:
##  $ ID   : chr  "JPN501" "JPN502" "JPN503" "JPN504" ...
##  $ Topic: chr  "sports" "education" "education" "sports" ...
##  $ Score: chr  "4" "4" "3" "4" ...
##  $ Type : int  134 0 0 160 0 0 121 0 0 139 ...
##  $ Token: int  638 0 0 712 0 0 402 0 0 520 ...
##  $ TTR  : num  0.42 NaN NaN 0.449 NaN ...
##  $ GI   : num  7.5 NaN NaN 8.48 NaN ...
##  $ MATTR: num  0.426 NaN NaN 0.455 NaN ...
##  $ AWL  : num  4.3 NaN NaN 4.23 NaN ...
##  $ ASL  : num  10.6 NaN NaN 12.3 NaN ...
##  $ NoS  : int  30 0 0 29 0 0 13 0 0 27 ...
##  $ Lang : num  2 2 2 2 2 2 2 2 2 2 ...
plot(all_indexes.df)


データフレームから、必要な情報だけを選び出す

特定の情報だけを集めたデータフレームを作る: data.frame(データフレーム\(見出し, データフレーム\)見出し, …)

例:IDとスコアと言語だけのデータセットを作る

ID_S_L <- data.frame(all_indexes.df$ID, all_indexes.df$Score, all_indexes.df$Lang)

head(ID_S_L)

条件に合うデータを取り出す:subset(データフレーム, 条件)

条件

  • 数値は同じか == 、大きいか > >= 、小さいか < <=
  • 文字は同じか == “文字列”" 、違うか != “文字列”
  • 「見出し」を使って表記

例: スコアが 6 のもの

subset(all_indexes.df, Score == 6)

例: スコアが 2以下 のもの

subset(all_indexes.df, Score <= 2)

複数の条件の組み合わせ

  • or 条件は |
  • and条件は &
  • ( ) を使って条件を入れ子にできる

例: スコアが 6 で、学習者のもの

subset(all_indexes.df, Score == 6 & Lang == 2)

例: スコアが 5 で、学習者のもの

subset(all_indexes.df, Score == 5 & Lang == 2)

例: スコアが 5 で、学習者のもので、トピックがsportsもしくはmoneyのもの

  • ( ) の使い方に注意
subset(all_indexes.df, Score == 5 & Lang == 2 & (Topic == "sports" | Topic == "money" ))
  • 別解
subset(all_indexes.df, Score == 5 & Lang == 2 & Topic != "education") 

条件に合うデータのうち、必要な情報だけのデータセットを作る: subset()のオプション select

  • 一つだけであれば select = 見出し
  • 複数ある場合は select = c(見出し1, 見出し2, …)

例:学習者データのスコアとTokenだけのデータセットを作る

L2token <- subset(all_indexes.df, Lang == 2, select = c(Score, Token)) 

str(L2token)
## 'data.frame':    1143 obs. of  2 variables:
##  $ Score: chr  "4" "4" "3" "4" ...
##  $ Token: int  638 0 0 712 0 0 402 0 0 520 ...
head(L2token)
plot(L2token)
## Warning in xy.coords(x, y, xlabel, ylabel, log): 強制変換により NA が生成されま
## した

plot(L2token$Token, L2token$Score)
## Warning in xy.coords(x, y, xlabel, ylabel, log): 強制変換により NA が生成されま
## した

型の変換: Scoreがchrになっている: as.numeric() で数値に変換

str(L2token)
## 'data.frame':    1143 obs. of  2 variables:
##  $ Score: chr  "4" "4" "3" "4" ...
##  $ Token: int  638 0 0 712 0 0 402 0 0 520 ...
L2token$Score <- as.numeric(L2token$Score)
## Warning: 強制変換により NA が生成されました
str(L2token)
## 'data.frame':    1143 obs. of  2 variables:
##  $ Score: num  4 4 3 4 4 3 4 3 4 3 ...
##  $ Token: int  638 0 0 712 0 0 402 0 0 520 ...
cor.test(L2token$Score, L2token$Token)
## 
##  Pearson's product-moment correlation
## 
## data:  L2token$Score and L2token$Token
## t = -0.57713, df = 1135, p-value = 0.564
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.07519020  0.04104952
## sample estimates:
##         cor 
## -0.01712822
plot(L2token$Token, L2token$Score)
abline(lm(L2token$Score~L2token$Token))

単回帰直線 Linear Model:要因が一つだけ

y = ax + b

lm(目的変数 ~ 説明変数, データ)

例:学習者データ jp_indexes.df

  1. 二種類のデータをプロットする (例えば、Score と Token)
  2. 単回帰直線のモデルを作る
  3. モデルをプロットする
  4. 概要を見る

注意: Scoreを数値型に変換しておく

jp_indexes.df$Score <- as.numeric(jp_indexes.df$Score)
## Warning: 強制変換により NA が生成されました

実行

plot(Score ~ Token, jp_indexes.df)
lm.result <- lm(Score ~ Token, jp_indexes.df)
abline(lm.result)

summary(lm.result)
## 
## Call:
## lm(formula = Score ~ Token, data = jp_indexes.df)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.5310 -0.5310 -0.4894  0.4690  1.5149 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  3.531e+00  2.708e-02 130.394   <2e-16 ***
## Token       -4.583e-05  7.942e-05  -0.577    0.564    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.7634 on 1135 degrees of freedom
##   ( 6 個の観測値が欠損のため削除されました )
## Multiple R-squared:  0.0002934,  Adjusted R-squared:  -0.0005874 
## F-statistic: 0.3331 on 1 and 1135 DF,  p-value: 0.564

重回帰分析 Multiple (Linear) Regression Analysis :要因が複数

http://sugiura-ken.org/wiki/wiki.cgi/exp?page=Multiple+Regression+Analysis

y = ax1 + bx2 + cx3 + ... + e

yが増える原因は複数ある。複数の要因が影響した結果がyとなる。

関数は lm()

lm(目的変数 ~ 説明変数1 + 説明変数2 + 説明変数3 + … , データ)

lm(Score ~ Token + Type + TTR + GI + MATTR + AWL + ASL + NoS, jp_indexes.df)
## 
## Call:
## lm(formula = Score ~ Token + Type + TTR + GI + MATTR + AWL + 
##     ASL + NoS, data = jp_indexes.df)
## 
## Coefficients:
## (Intercept)        Token         Type          TTR           GI        MATTR  
##    3.754562    -0.001634     0.026017    -7.406015    -0.585338    10.314709  
##         AWL          ASL          NoS  
##    0.052445     0.006108     0.003172

注意点1:目的変数が正規分布しているか。

重回帰分析の注意点

  • データが前提としている確率分布(この場合は正規分布)に従っているか

  • 重回帰分析では、目的変数の誤差(説明変数からのズレ)が正規分布に従っているという前提

  • Normal Q-Q プロットで確認 ** Y軸:観測データの標準化残差 ** X軸:期待値

hist(jp_indexes.df$Score)

qqnorm(jp_indexes.df$Score)   # 散布図
qqline(jp_indexes.df$Score)   # 直線

注意点2:変数の規模・単位が違うので、影響を等しく比べられるように「標準化」する: scale()

lm(Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + scale(MATTR) + scale(AWL) + scale(ASL) + scale(NoS), jp_indexes.df)
## 
## Call:
## lm(formula = Score ~ scale(Token) + scale(Type) + scale(TTR) + 
##     scale(GI) + scale(MATTR) + scale(AWL) + scale(ASL) + scale(NoS), 
##     data = jp_indexes.df)
## 
## Coefficients:
##  (Intercept)  scale(Token)   scale(Type)    scale(TTR)     scale(GI)  
##      1.85286      -0.46540       1.62257      -0.48619      -0.56766  
## scale(MATTR)    scale(AWL)    scale(ASL)    scale(NoS)  
##      0.66889       0.01772       0.01689       0.03540

注意点3: 説明する要因間に相関が高いものがあるので、 VIF(Variance Inflation factor)をチェック

  • それらの要因の元になっている潜在的な要因がある可能性がある
  • 複数の相関が高いものを重ねて入れてしまうと、その影響が重複して評価されてしまう。
  • 10以上はダメ、2以下がよい
  • 組み合わせによるので、試行錯誤が必要。
  • carパッケージをインストールし、library(car)をする。
  • 関数はvif()

パッケージのインストール: install.packages(“パッケージ名”)

パッケージを使えるようにする: library(パッケージ名)

  • ちょと時間がかかる。
install.packages("car", repos="http://cran.rstudio.com/")
library(car)
##  要求されたパッケージ carData をロード中です

標準化したデータで重回帰分析をした結果を保存

jukaiki.result <- lm(Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + scale(MATTR) + scale(AWL) + scale(ASL) + scale(NoS), jp_indexes.df)

重回帰した結果をVIFでチェック

vif(jukaiki.result)
## scale(Token)  scale(Type)   scale(TTR)    scale(GI) scale(MATTR)   scale(AWL) 
##   100.530629   290.802743   479.778827   135.493414   517.493991     1.127758 
##   scale(ASL)   scale(NoS) 
##    11.980648    25.996185
  • 10以上はダメ、2以下がよい
  • 重複していると考えられるものを除いて、モデルと作って試す。

語彙の多様性指標はMATTRに限定

jukaiki.result2 <- lm(Score ~ scale(Token) + scale(Type) + scale(MATTR) + scale(AWL) + scale(ASL) + scale(NoS), jp_indexes.df)
vif(jukaiki.result2)
## scale(Token)  scale(Type) scale(MATTR)   scale(AWL)   scale(ASL)   scale(NoS) 
##    41.450953    20.382691     6.918761     1.119840    11.192397    24.187514

Typeを外してみる

jukaiki.result3 <- lm(Score ~ scale(Token) + scale(MATTR) + scale(AWL) + scale(ASL) + scale(NoS), jp_indexes.df)
vif(jukaiki.result3)
## scale(Token) scale(MATTR)   scale(AWL)   scale(ASL)   scale(NoS) 
##    24.484244     1.764875     1.119781    10.524400    22.445986

文の数を外してみる

jukaiki.result4 <- lm(Score ~ scale(Token) + scale(MATTR) + scale(AWL) + scale(ASL), jp_indexes.df)
vif(jukaiki.result4)
## scale(Token) scale(MATTR)   scale(AWL)   scale(ASL) 
##     1.815905     1.757304     1.118768     1.207492
  • すべて2以下

得られたモデル

lm(Score ~ scale(Token) + scale(MATTR) + scale(AWL) + scale(ASL), jp_indexes.df)
## 
## Call:
## lm(formula = Score ~ scale(Token) + scale(MATTR) + scale(AWL) + 
##     scale(ASL), data = jp_indexes.df)
## 
## Coefficients:
##  (Intercept)  scale(Token)  scale(MATTR)    scale(AWL)    scale(ASL)  
##     3.317656      0.124462      0.025752      0.017175      0.003896
  • 係数が大きいものほど影響力が強い
    • TokenとAWL

結果を詳しく見る:summary()

summary(jukaiki.result4)
## 
## Call:
## lm(formula = Score ~ scale(Token) + scale(MATTR) + scale(AWL) + 
##     scale(ASL), data = jp_indexes.df)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.5874 -0.4644 -0.3504  0.5251  1.6669 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  3.317656   0.098502  33.681   <2e-16 ***
## scale(Token) 0.124462   0.070579   1.763   0.0786 .  
## scale(MATTR) 0.025752   0.045873   0.561   0.5749    
## scale(AWL)   0.017175   0.036602   0.469   0.6392    
## scale(ASL)   0.003896   0.038026   0.102   0.9184    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6746 on 376 degrees of freedom
##   ( 762 個の観測値が欠損のため削除されました )
## Multiple R-squared:  0.01157,    Adjusted R-squared:  0.001058 
## F-statistic: 1.101 on 4 and 376 DF,  p-value: 0.356
  • TokenとAWLに加え、ASLも有意
  • MATTRはなくてもよさそう
  • 説明率R^2 (Adjusted R-squared: 0.7304)が、このモデルの説明率を示している。

MATTRが重要でないなら、外して、モデルを作り直してみる

jukaiki.result5 <- lm(Score ~ scale(Token) + scale(AWL) + scale(ASL), jp_indexes.df)
summary(jukaiki.result5)
## 
## Call:
## lm(formula = Score ~ scale(Token) + scale(AWL) + scale(ASL), 
##     data = jp_indexes.df)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.5892 -0.4626 -0.3412  0.5205  1.6261 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   3.34800    0.08227  40.697   <2e-16 ***
## scale(Token)  0.10124    0.05713   1.772   0.0772 .  
## scale(AWL)    0.02244    0.03535   0.635   0.5258    
## scale(ASL)    0.00315    0.03797   0.083   0.9339    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.674 on 377 degrees of freedom
##   ( 762 個の観測値が欠損のため削除されました )
## Multiple R-squared:  0.01074,    Adjusted R-squared:  0.002872 
## F-statistic: 1.365 on 3 and 377 DF,  p-value: 0.2532
  • 説明率が 0.7311となった。むしろこの方が高い。

説明変数に何を入れるかを自動で行う変数選択:step()

  • 説明変数をすべて入れたモデル

  • 欠損値があるとsetp関数がエラーになるので、na.omit()でデータフレームを処置しておく

jukaiki.result <- step(lm(Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + scale(MATTR) + scale(AWL) + scale(ASL) + scale(NoS), na.omit(jp_indexes.df)))
## Start:  AIC=-289.3
## Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + 
##     scale(MATTR) + scale(AWL) + scale(ASL) + scale(NoS)
## 
##                Df Sum of Sq    RSS     AIC
## - scale(NoS)    1   0.00718 170.08 -291.28
## - scale(ASL)    1   0.00905 170.09 -291.28
## - scale(AWL)    1   0.10578 170.18 -291.06
## - scale(TTR)    1   0.18722 170.26 -290.88
## - scale(MATTR)  1   0.32854 170.41 -290.56
## - scale(Token)  1   0.35740 170.43 -290.50
## <none>                      170.08 -289.30
## - scale(GI)     1   0.90373 170.98 -289.28
## - scale(Type)   1   0.92183 171.00 -289.24
## 
## Step:  AIC=-291.28
## Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + 
##     scale(MATTR) + scale(AWL) + scale(ASL)
## 
##                Df Sum of Sq    RSS     AIC
## - scale(ASL)    1   0.00216 170.09 -293.27
## - scale(AWL)    1   0.10639 170.19 -293.04
## - scale(TTR)    1   0.19912 170.28 -292.83
## - scale(MATTR)  1   0.33451 170.42 -292.53
## - scale(Token)  1   0.41521 170.50 -292.35
## <none>                      170.08 -291.28
## - scale(GI)     1   0.91156 171.00 -291.24
## - scale(Type)   1   0.91703 171.00 -291.23
## 
## Step:  AIC=-293.27
## Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + 
##     scale(MATTR) + scale(AWL)
## 
##                Df Sum of Sq    RSS     AIC
## - scale(AWL)    1   0.11612 170.20 -295.01
## - scale(TTR)    1   0.20182 170.29 -294.82
## - scale(MATTR)  1   0.33825 170.42 -294.52
## - scale(Token)  1   0.41319 170.50 -294.35
## <none>                      170.09 -293.27
## - scale(GI)     1   0.91286 171.00 -293.24
## - scale(Type)   1   0.91745 171.00 -293.23
## 
## Step:  AIC=-295.01
## Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + 
##     scale(MATTR)
## 
##                Df Sum of Sq    RSS     AIC
## - scale(TTR)    1   0.23031 170.43 -296.50
## - scale(MATTR)  1   0.37857 170.58 -296.17
## - scale(Token)  1   0.39513 170.60 -296.13
## - scale(GI)     1   0.89184 171.09 -295.02
## <none>                      170.20 -295.01
## - scale(Type)   1   0.89587 171.10 -295.01
## 
## Step:  AIC=-296.5
## Score ~ scale(Token) + scale(Type) + scale(GI) + scale(MATTR)
## 
##                Df Sum of Sq    RSS     AIC
## - scale(Token)  1   0.32542 170.76 -297.77
## - scale(MATTR)  1   0.59637 171.03 -297.17
## - scale(GI)     1   0.76369 171.20 -296.80
## - scale(Type)   1   0.77643 171.21 -296.77
## <none>                      170.43 -296.50
## 
## Step:  AIC=-297.77
## Score ~ scale(Type) + scale(GI) + scale(MATTR)
## 
##                Df Sum of Sq    RSS     AIC
## - scale(MATTR)  1   0.40796 171.17 -298.86
## - scale(GI)     1   0.46730 171.22 -298.73
## - scale(Type)   1   0.75596 171.51 -298.09
## <none>                      170.76 -297.77
## 
## Step:  AIC=-298.86
## Score ~ scale(Type) + scale(GI)
## 
##               Df Sum of Sq    RSS     AIC
## - scale(GI)    1   0.06679 171.23 -300.71
## <none>                     171.17 -298.86
## - scale(Type)  1   0.94503 172.11 -298.77
## 
## Step:  AIC=-300.72
## Score ~ scale(Type)
## 
##               Df Sum of Sq    RSS     AIC
## <none>                     171.23 -300.71
## - scale(Type)  1    1.8699 173.10 -298.58

AICが一番小さいものを最善のモデルとして選択

summary(jukaiki.result)
## 
## Call:
## lm(formula = Score ~ scale(Type), data = na.omit(jp_indexes.df))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.5755 -0.4690 -0.3582  0.5201  1.6657 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  3.48031    0.03444 101.067   <2e-16 ***
## scale(Type)  0.07015    0.03448   2.034   0.0426 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6722 on 379 degrees of freedom
## Multiple R-squared:  0.0108, Adjusted R-squared:  0.008192 
## F-statistic: 4.139 on 1 and 379 DF,  p-value: 0.04261
  • これで説明率 0.76 まで上がった。
  • 説明変数は、TTR, GI, AWL, ASL, NoS
jukaiki.best.result <- lm(formula = Score ~ scale(TTR) + scale(GI) + scale(AWL) + scale(ASL) + 
    scale(NoS), data = jp_indexes.df)

vif(jukaiki.best.result)
## scale(TTR)  scale(GI) scale(AWL) scale(ASL) scale(NoS) 
##  10.076907   8.136216   1.116245   5.460248  11.455641
  • VIFは10をこえないこと、というのが原則。これくらいならよしとするか。

回帰診断図(分析の評価)

par(mfrow=c(2,2),oma = c(1,1,2,1),mar = c(4, 4, 2, 1))
plot(jukaiki.best.result)

  1. Residuals vs Fitted
    • 残差と予測値:予測値によって変な偏りがないか:赤い横線が水平のほうがよい
  2. Normal Q-Q
    • 残差が正規分布に従っている方が良い:線の上にまっすくに乗っているか?
  3. Scale-Location
    • 標準化残差の平方根と予測値:残差を標準化したもの(Y軸):赤い横線が水平のほうがよい
  4. Residuals vs Leverage
    • 標準化残差とてこ比:影響力(Cookの距離)の大きさのばらつきを見る:0.5の中に納まっている方が良い