出力（プログラムで処理した結果）をファイルに保存する

コンソール（Rの画面上）に表示される結果を、ファイルに保存する。
出力先を変えることを「リダイレクト」と呼ぶ。

sink()

sink(“ファイル名”) とすることで、出力がファイルに保存される。処理の前に書いておく。
sink() と書くことで、出力先がコンソールになる。処理の後に書いておく。
- これを書いておかないと、それ以降、Rの処理結果が画面に表示されなくなる。

例：myVoc() と myVocSink()

myVoc()

myVoc <- function(){

    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))
        
        cat(length(type), length(token), ttr, gi, "\n")
    }

}

myVocSink()

出力先のファイル名の前に、上の階層を指定する（ ../ を前につける）点に注意
myVocのプログラムが、作業ディレクトリー内のすべてのファイルを処理するものなので、
出力結果のファイルは、その中に入れずに、外（上の階層）に出しておくようにする。

myVocSink <- function(){
    sink("../myVocResult.txt")
    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))
        
        cat(length(type), length(token), ttr, gi, "\n")
    }
    sink()
}

出力（プログラムで処理した結果）を「データフレーム」に保存する

データの種類

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

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

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

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

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

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

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

ベクトルの入れ物を用意（初期化）
それぞれのベクトルに数値を入れていく
ベクトルをまとめてデータフレームにする

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_2/NICER_NS")

myIndexesDf()

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

setwd("../NICER1_3_2/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)

♪　演習

演習１：TokenとTTRの関係を散布図に

演習２：TokenとGIの関係を散布図に

解答例１

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()

♪　演習

演習３：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, 追加するベクトルデータ)

cbind(もとのdf, 追加するベクトルデータ) でもできる

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

setwd("../NICER1_3_2/NICER_NS")

mattr <- mattr2()

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

head(indexes.df2)

演習４：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のスコア)を取り込んでみる

考え方

検索（lines.tmpに一文一行でデータが入っているとして）

criterion.tmp <- grep("@Criterion", lines.tmp, value = T)

置換・削除：その行の値のところだけ残して前を消す

score <- gsub("@Criterion:\t", "", criterion.tmp)

♪　演習

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

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

(^^♪ 解答例

setwd("../NICER1_3_2/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()

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

setwd("../NICER1_3_2/NICER_NS")

score <- myScore()


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

head(indexes.df3)

演習７：参加者ID（JPN501とかNS501）も、データフレームに追加する

setwd("../NICER1_3_2/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_2/NICER_NS")

id <- myID()


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

head(indexes.df4)

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

plot(indexes.df4)

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

それができれば、学習者データのデータフレームもすぐに出ます。

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_2/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_2/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(データフレーム１, データフレーム２)

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':    452 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 160 121 139 175 124 151 98 104 99 ...
##  $ Token: int  638 712 402 520 840 522 724 396 526 366 ...
##  $ TTR  : num  0.42 0.449 0.602 0.535 0.417 ...
##  $ GI   : num  7.5 8.48 8.53 8.62 8.54 ...
##  $ MATTR: num  0.426 0.455 0.606 0.539 0.423 ...
##  $ AWL  : num  4.3 4.23 4.75 4.77 4 ...
##  $ ASL  : num  10.63 12.28 15.46 9.63 16.8 ...
##  $ NoS  : int  30 29 13 27 25 20 26 20 19 14 ...
##  $ 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(見出し１, 見出し２, …)

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

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

str(L2token)

## 'data.frame':    381 obs. of  2 variables:
##  $ Score: chr  "4" "4" "3" "4" ...
##  $ Token: int  638 712 402 520 840 522 724 396 526 366 ...

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':    381 obs. of  2 variables:
##  $ Score: chr  "4" "4" "3" "4" ...
##  $ Token: int  638 712 402 520 840 522 724 396 526 366 ...

L2token$Score <- as.numeric(L2token$Score)

## Warning: 強制変換により NA が生成されました

str(L2token)

## 'data.frame':    381 obs. of  2 variables:
##  $ Score: num  4 4 3 4 4 3 4 3 4 3 ...
##  $ Token: int  638 712 402 520 840 522 724 396 526 366 ...

cor.test(L2token$Score, L2token$Token)

## 
##  Pearson's product-moment correlation
## 
## data:  L2token$Score and L2token$Token
## t = 27.644, df = 377, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.7820537 0.8490645
## sample estimates:
##       cor 
## 0.8183211

plot(L2token$Token, L2token$Score)
abline(lm(L2token$Score~L2token$Token))

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

y = ax + b

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

因果関係が想定されている。xが増えると、その結果yが増える。
yが増える要因（原因）は一つxだけ。（単回帰分析）

例：学習者データ jp_indexes.df

二種類のデータをプロットする　（例えば、Score と Token）
単回帰直線のモデルを作る
モデルをプロットする
概要を見る

注意： 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 
## -1.5569 -0.2748  0.0109  0.2933  1.3232 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1.6669449  0.0708166   23.54   <2e-16 ***
## Token       0.0033220  0.0001202   27.64   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4396 on 377 degrees of freedom
##   ( 2 個の観測値が欠損のため削除されました )
## Multiple R-squared:  0.6696, Adjusted R-squared:  0.6688 
## F-statistic: 764.2 on 1 and 377 DF,  p-value: < 2.2e-16

重回帰分析 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(目的変数 ~ 説明変数１ + 説明変数２ + 説明変数３ + … , データ)

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  
##   -0.535503     0.001256    -0.008998    -9.016308     0.577063     2.682524  
##         AWL          ASL          NoS  
##    0.506956     0.041955     0.013926

係数の大きさが、影響の大きさと見なされる

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

重回帰分析の注意点

データが前提としている確率分布（この場合は正規分布）に従っているか
重回帰分析では、目的変数の誤差（説明変数からのズレ）が正規分布に従っているという前提
Normal Q-Q プロットで確認 ** Y軸：観測データの標準化残差 ** X軸：期待値

hist(jp_indexes.df$Score)

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

注意点２：変数の規模・単位が違うので、影響を等しく比べられるように「標準化」する： 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)  
##      3.52043       0.23632      -0.29043      -0.59159       0.55935  
## scale(MATTR)    scale(AWL)    scale(ASL)    scale(NoS)  
##      0.17387       0.17130       0.11596       0.09719

注意点３：説明する要因間に相関が高いものがあるので、 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) 
##   101.018056   291.500537   474.789406   135.912033   512.414581     1.125987 
##   scale(ASL)   scale(NoS) 
##    12.071882    25.998623

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.462335    20.388950     6.864829     1.117988    11.244092    24.131301

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.554126     1.756094     1.117927    10.574447    22.395073

文の数を外してみる

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.812922     1.750001     1.117058     1.210238

すべて２以下

得られたモデル

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.520855      0.625732      0.006836      0.176912      0.051028

係数が大きいものほど影響力が強い
- 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 
## -1.38673 -0.25294 -0.00528  0.25302  1.15642 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  3.520855   0.020380 172.763  < 2e-16 ***
## scale(Token) 0.625732   0.027469  22.779  < 2e-16 ***
## scale(MATTR) 0.006836   0.027054   0.253   0.8006    
## scale(AWL)   0.176912   0.021538   8.214 3.53e-15 ***
## scale(ASL)   0.051028   0.022397   2.278   0.0233 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3967 on 374 degrees of freedom
##   ( 2 個の観測値が欠損のため削除されました )
## Multiple R-squared:  0.7331, Adjusted R-squared:  0.7302 
## F-statistic: 256.8 on 4 and 374 DF,  p-value: < 2.2e-16

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 
## -1.37835 -0.25224 -0.00619  0.25823  1.16749 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   3.52086    0.02035 172.979  < 2e-16 ***
## scale(Token)  0.62168    0.02227  27.911  < 2e-16 ***
## scale(AWL)    0.17830    0.02080   8.571 2.72e-16 ***
## scale(ASL)    0.05083    0.02236   2.274   0.0236 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3963 on 375 degrees of freedom
##   ( 2 個の観測値が欠損のため削除されました )
## Multiple R-squared:  0.733,  Adjusted R-squared:  0.7309 
## F-statistic: 343.2 on 3 and 375 DF,  p-value: < 2.2e-16

説明率が 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=-740.94
## Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + 
##     scale(MATTR) + scale(AWL) + scale(ASL) + scale(NoS)
## 
##                Df Sum of Sq    RSS     AIC
## - scale(MATTR)  1    0.0222 51.187 -742.77
## - scale(Type)   1    0.1096 51.275 -742.13
## - scale(NoS)    1    0.1373 51.303 -741.92
## - scale(Token)  1    0.2091 51.374 -741.39
## <none>                      51.165 -740.94
## - scale(TTR)    1    0.2774 51.443 -740.89
## - scale(ASL)    1    0.4230 51.588 -739.82
## - scale(GI)     1    0.8721 52.037 -736.53
## - scale(AWL)    1    9.8783 61.044 -676.03
## 
## Step:  AIC=-742.77
## Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + 
##     scale(AWL) + scale(ASL) + scale(NoS)
## 
##                Df Sum of Sq    RSS     AIC
## - scale(NoS)    1    0.1434 51.331 -743.71
## - scale(Type)   1    0.1467 51.334 -743.69
## - scale(Token)  1    0.2294 51.417 -743.08
## <none>                      51.187 -742.77
## - scale(ASL)    1    0.4373 51.625 -741.55
## - scale(GI)     1    1.0482 52.236 -737.09
## - scale(TTR)    1    3.7103 54.898 -718.25
## - scale(AWL)    1   10.0173 61.205 -677.04
## 
## Step:  AIC=-743.71
## Score ~ scale(Token) + scale(Type) + scale(TTR) + scale(GI) + 
##     scale(AWL) + scale(ASL)
## 
##                Df Sum of Sq    RSS     AIC
## - scale(Type)   1    0.2030 51.534 -744.22
## <none>                      51.331 -743.71
## - scale(Token)  1    0.6320 51.963 -741.08
## - scale(ASL)    1    0.9065 52.237 -739.08
## - scale(GI)     1    1.3288 52.660 -736.03
## - scale(TTR)    1    4.8778 56.209 -711.31
## - scale(AWL)    1   10.0497 61.380 -677.95
## 
## Step:  AIC=-744.22
## Score ~ scale(Token) + scale(TTR) + scale(GI) + scale(AWL) + 
##     scale(ASL)
## 
##                Df Sum of Sq    RSS     AIC
## <none>                      51.534 -744.22
## - scale(Token)  1    0.8577 52.392 -739.96
## - scale(ASL)    1    0.9203 52.454 -739.51
## - scale(TTR)    1    5.9727 57.507 -704.66
## - scale(GI)     1    7.3450 58.879 -695.72
## - scale(AWL)    1   10.2535 61.787 -677.44

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

summary(jukaiki.result)

## 
## Call:
## lm(formula = Score ~ scale(Token) + scale(TTR) + scale(GI) + 
##     scale(AWL) + scale(ASL), data = na.omit(jp_indexes.df))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.32107 -0.24903 -0.01515  0.25641  1.49021 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   3.52243    0.01909 184.489  < 2e-16 ***
## scale(Token)  0.16837    0.06758   2.492   0.0132 *  
## scale(TTR)   -0.39965    0.06078  -6.575 1.64e-10 ***
## scale(GI)     0.39869    0.05468   7.291 1.85e-12 ***
## scale(AWL)    0.17389    0.02019   8.615  < 2e-16 ***
## scale(ASL)    0.05430    0.02104   2.581   0.0102 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3717 on 373 degrees of freedom
## Multiple R-squared:  0.7663, Adjusted R-squared:  0.7632 
## F-statistic: 244.7 on 5 and 373 DF,  p-value: < 2.2e-16

これで説明率 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.013866   8.159103   1.114485   5.472266  11.398238

VIFは10をこえないこと、というのが原則。これくらいならよしとするか。

回帰診断図（分析の評価）

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

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

参考：「Rと回帰分析」（https://www1.doshisha.ac.jp/~mjin/R/Chap_14/14.html）

R for Learner Corpus Research 2021: 05

sugiura

2021/12/03

出力（プログラムで処理した結果）をファイルに保存する

sink()

例：myVoc() と myVocSink()

myVoc()

myVocSink()

出力（プログラムで処理した結果）を「データフレーム」に保存する

データの種類

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

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

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

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

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

実行

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

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

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

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

データフレームの一部（特定の列（カラム））だけ取り出す： データフレーム$見出し

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

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

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

基本的なグラフ

boxplot()

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

散布図： plot()

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

棒グラフ： barplot()

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

ヒストグラム： hist()

♪ 演習

演習１：TokenとTTRの関係を散布図に

演習２：TokenとGIの関係を散布図に

解答例１

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

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

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

TTR, GIの問題点と解決案

Covington and McFall (2008, 2010)

MATTR

♪ 演習

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

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

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

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

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

ヘッダー情報の取り込み

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

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

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

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

考え方

♪ 演習

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

(^^♪ 解答例

実行

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

演習７：参加者ID（JPN501とかNS501）も、データフレームに追加する

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

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

実行

母語話者データ

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

学習者データ

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

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

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

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

データ構造の確認

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

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

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

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

条件

例： スコアが 6 のもの

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

複数の条件の組み合わせ

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

データフレームの一部（特定の列（カラム））だけ取り出す：データフレーム$見出し

♪　演習

♪　演習

♪　演習

例：スコアが 6 のもの

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

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

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

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

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

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