EcoPlate解析 (1)
１）エコプレートのデータファイルの準備
マイクロプレートリーダーにはそのモデルごとに読み取った生データを変換するためのアプリケーションが付属している場合が多いですが、そのような前処理は全く必要なく、以下のような生データをテキストファイルかもしくはcsvファイルで自分のコンピュータに保存します。
csvファイルなら、エクセルで開くと以下の感じです。
どちらの形式でも構いませんが、以下はcsvファイルを使うこととします。このようなファイルを、ひとつの研究プロジェクトでは、複数のサンプルごと、ひとつのサンプルの培養時間ごとに数多く作ることになります。ファイルの準備で重要なことは以下の2つです。
１）フォルダ名・ファイル名に一貫性を持たせること
たとえば、"サンプル日時-処理名orサンプル地点-測定時間.csv"などのような感じです。ファイル名は数字とアルファベットだけ使いましょう。またスペースをファイル名の中に含めるのはやめましょう。
２）ファイルの中身のフォーマットに一貫性を持たせること
使用するマイクロプレートリーダーごとに生成されるファイルのフォーマットは異なります。通常、生データの上に測定日時や測定プロトコルなどのメタデータが含まれますが（たとえば、上のcsvファイルの１から6行目）、それは手作業ではいじらずにそのままにします。また空欄の行などもそのままにします。これにより同一のマイクロプレートリーダーを使っている限り生データのフォーマットに一貫性が保たれます。
２）Rスクリプト全体のデザイン
Rに限らず言語によらずスクリプトやソースコードの可読性を高く維持することは最も重要なことであると同時に実現するのが難しいことでもあります。基本的には以下のように、スクリプトの種類別に「ブロック」を作ってデザインするのがよいと思います。
##########スクリプト全体へのコメント###################必要なライブラリを全部読み込む（ライブラリ・ブロック）##########library(library_name1)library(library_name2)...########必要な関数をすべて定義する（関数ブロック）#############func1<-function(...){}func2<-function(...){}...########関数などを読み出して計算を実行する（解析ブロック）############解析の種類ごとにブロックを作る##################Analysis01#######a<-func1(...)####Analysis02#######....#######関数をつくるまでにやった予備計算のスクリプトなどを最後に残しておいてもよい####おわりしかし、今回の講習会で配布するサンプルコードは解析方法を学習するためのものですので、上のような構造にはなっていません。サンプルコードは、以下のような順序で書かれています。
#1)解析1に必要なライブラリを読み込む（もしあれば）#2)解析に必要な関数を定義する前に予備的計算を行う#3)関数を定義する#4)定義した関数を使って解析を行う#1)解析2に必要なライブラリを読み込む（以下2)以降に戻って、解析の種類だけ繰り返す）講習後に自分でスクリプトを書くときは、構造に気をつけて書くのがよいと思います。
３）エコプレートのデータファイルの読み込み
"サンプル日時-処理名orサンプル地点-測定時間.csv"のように一貫性をもってデータファイルに名前が付いているとします。この一貫性があれば、複数のデータファイルを自動的にすべて読み込むようなスクリプトを書くことが可能です。そのようなスクリプトを関数として定義するための準備として、複数の文字列をくっつけて新しい文字列を作るための関数"paste"の使い方に慣れる必要があります。
以下の例（の最終行）は、ACGCTTT, AGGTTC, data01-021,という3つの文字列を、間に"_"を挟んで（sep="_"というオプションで）くっつけて、ACGCTT_AGCCTTC_data01-023という文字列を生成するためのスクリプトです。
####Preparation#########use the function "paste" to generate the sequence of relative path to efficiently access to multiple files of ecoplate data#example 01paste("ACGCTT", "AGCCTTC", "data01-023", sep="_")paste()関数の中に、直接長い文字列を書いていくと可読性が下がるので、以下のように複数の行に分けて書くのもよいです。
place_name<-"osaka"file_name<-"sample01.csv"file_path<-paste("./data", place_name, file_names, sep="/")file_pathこの例を応用して、実際のファイルをいくつか読み込むスクリプトを書いてみます。ここでは、親フォルダの中にフォルダが二つあり(data_merged, r_analysis)、Rスクリプトはr_analysisの中にあり、データファイルはdata_mergedの中に保存されているとします。データファイルは複数あって、それらの名前は、20150512-1500higashihori_24.csv, 20150512-1500higashihori_48.csv, 20150512-1500higashihori_72.csv, としましょう。フォーマットとしては、"サンプル採種年月日-サンプリング時刻サンプリング場所_エコプレート培養時間.csv"となっています。
この場合、r_analysisの中からひとつ上の階層に出てから(../)、data_mergedにアクセスするように、データファイルへの相対パスを書く必要があります。つまり、
../data_merged/20150512-1500higashihori_XX.csv
という感じです。XXは24、48、72のどれかです。XXの手前までの情報として以下の2行をスクリプトに書きます。
############Simple trial to load three files at the same time and stock the data into list##############path_file<-"../data_merged/" #folder namedata_name<-"20150512-1500higashihori" #data name培養時間XXの部分は順次変えないといけないので、ベクトルとして用意します。
time_m<-c(24,48,72) #sequence of time_pointここで作ったtime_mというオブジェクトの各要素へのアクセスとデータのタイプの確認は以下のようにできます。
#check the typetime_m[1]#check the typeclass(time_m[1]) #--> this should be converted into character to be used in the part of file pathtime_m[1]のクラスは数値(numeric)なので、ファイルへの相対パスの一部として文字列に組み込むには文字列に変換する関数as.characterを使う必要があります。次のスクリプトは、まずファイルをひとつ（20150512-1500higashihori_24.csv ）だけ読み込むための相対パスの生成のためのものです。
#To generate file path "../higashi_hori_bashi_rain/20150512-1500higashihori_24.csv""j<-1file_ph<-paste(path_file, data_name, "_", as.character(time_m[j]), ".csv", sep="")file_phfile_phが"../data_merged/20150512-1500higashihori_24.csv"となっているのが分かると思います。
それでは次にまず、R入門のページで説明した「リスト」を、複数のエコプレートデータを読み込むためのハコとして用意します。
#generate empy list ecoplate data are saved as a listdata_box<-list()つづいて、このdata_boxの一個目の要素data_box[[1]]に20150512-1500higashihori_24.csvの内容を読み込むことを考えます。このファイルの中身は上のエクセルの図にあるように、上から6行目まではメタデータ、7行目はエコプレート本体のウェルの通し番号（１，２、．．．１２）がついていてこれらはすべて読み込む必要がありません。また、一列目には同じくメタデータの一部と、エコプレート本体のウェルの通し番号（A-H）がついていてこれは読み込まざるを得ませんが読み込んだ後に削除します。最初の７行を読み込むときにスキップするとして、以下のようにスキップ数を指定します。
#the first several lines in csv file include metadata, but should be skippedno_skip<-77行分を飛ばすといきなりエコプレートのデータが始まるので、ヘッダーも必要なく(header=FALSE)、データを読み込み、中身を見てみましょう。
#read the text file and saved into the item of the listdata_box[[j]] <- read.csv(file_ph, skip=no_skip, header=FALSE) View(data_box[[j]])View()関数で中身を見てみると確かに、列名(colname)にはV1-V13と意味のない記号が割り当てられています。V1の列はA-Hが記録されているだけで必要ないので削除し、かつ同じオブジェクト(data_box[[j])に上書きしましょう。data_box[[j]][, -1]の"-1"の部分が、1列目を取り除くという意味になります。
#Deleve the first column (A,B,...H)data_box[[j]] <-data_box[[j]][,-1] つづいて分かりやすいように行（row）と列（column）に名前をつけましょう。ここでは、エコプレート本体のウェルとの対応が分かればよいので簡単な名前で十分です。物質名を入れる必要はありません。
#set namesrownames(data_box[[j]])<-c("A", "B", "C", "D", "E", "F", "H", "G") colnames(data_box[[j]])<-c("V1","V2","V3", "V4","V5","V6","V7","V8", "V9", "V10", "V11", "V12") これで、View(data_box[[j]])を実行すれば順序良くデータが読み込まれたことが確認できるはずです（上のエクセルのキャプチャー図とは値が一致しなくてよいので注意）。
以上で、ひとつのデータファイルの中身ををリストに読み込む過程が理解できたはずです。今度はこれをループを使って自動化していきましょう。time_mには、time_m<-c(24,48,72) によりすでに3つの値が読み込まれ3つのファイルからのデータを読み込む準備ができています。したがってforループを以下のように使えば、3つのファイルの中身を順番にリストの要素各々へ読み込むことができます。forループの中身はこれまで説明したスクリプトが組み合わされています。
#If you like to load three files ***24.csv, ***48.csv, ***72.csv at the same time, you can use the loop functionfor(j in 1:3) { #file path file_ph<-paste(path_file, data_name, "_", as.character(time_m[j]), ".csv", sep="") data_box[[j]] <- read.csv(file_ph, skip=no_skip, header=FALSE) #Deleve the first column (A,B,...H) data_box[[j]] <-data_box[[j]][,-1] #set names rownames(data_box[[j]])<-c("A", "B", "C", "D", "E", "F", "G", "H") colnames(data_box[[j]])<-c("V1","V2","V3", "V4","V5","V6","V7","V8", "V9", "V10", "V11", "V12") }うまくいったかどうかは、View(data_box[[2]]) などで中身を確かめればチェックできます。
つぎに、このforループをコアとした関数を用意し、より多くのファイルをファイル名などを指定して読み込める形にします。関数を定義する前に、以下のようなエラー処理について考えます。time_m2というベクトルを用意し、実際には存在しない培養時間（28）のファイルにもアクセスしようとするとどうなるか試してみます。
#Add one more complexity to avoid the error#Skip the index j if the file does not existtime_m2<-c(24,28, 48,72) #sequence of time_point, noting that the data (file) at 28 does not existtime_m2の一つ目の要素、time_m2[1]を使えば、確かに存在するファイルにアクセスを試みるので、以下のスクリプトで"e"というオブジェクトのクラスはは、確かにデータを読み込んだ"data.frame"になります。データも確かに読み込まれています。
#check what happens without or with error by class(e)#without errorj<-1file_ph<-paste(path_file, data_name, "_", as.character(time_m2[j]), ".csv", sep="")e<-try(read.csv(file_ph, skip=no_skip, header=T), silent=FALSE) #error managementclass(e) eこんどは、time_m2[2]]を使って、存在しないファイルにアクセスしたときの結果を見てみます。
#with errorj<-2file_ph<-paste(path_file, data_name, "_", as.character(time_m2[j]), ".csv", sep="")e<-try(read.csv(file_ph, skip=no_skip, header=T), silent=FALSE) #error managementclass(e) これらの行を実行すると、コンソールにエラーメッセージが表示されるとともにeのクラスは"try-error"となります。
さて、このtry()という関数の中に実際実行したいことを書いて（ここではread.csv）関数の中に組み込むと、たとえエラーが出てもそこで計算が止まらずエラーが出ないスクリプトに進むことができるようになります。いきなり複雑な見た目になりますが、以上のことを踏まえると、以下のような関数が定義できます。最初の数行は関数の定義に関する情報です。エラー処理の部分だけ強調のため太字にしてあります。if(class(e)=="try-error") next というエラーチェックの条件文によって、ファイルが存在せず一行上のread.csv()がエラーを返したときは、else以降で実際にファイルを読んだり名前を変えたりするプロセスを飛ばして、次のループインデックスiに進むという構造になっています。
#funtion version 1 using file name###############list of parameters########################relative_path & folder_name: need to specify the folder position of data #date_time_sample: the string to specify the sample date and time, which should have the correspondance with thte file name#place_sample: the name of sample place, which is the part of file name#measured_time: the vector of the measured time#no_skip_line: The first some lines of the raw output from microplate reader represent metadata, which should be excluded##################Output is a list; each element of the list is the data frame to have ecoplate 96 valuesload_ecoplate_data_osaka <- function(relative_path="../", folder_name="data_merged", date_time_sample="20150115-0700", place_sample ="higashihori", measured_time=c(24,48,72), no_skip_line=7) #The parameters that already have values in the functional definition are recognized as default setting.{ data_list <-list() #generate empty list, ecoplate data are saved as a list loop_length=length(measured_time) #count the length of loop by the number of files for(i in 1:loop_length) { file_name <- paste(relative_path,folder_name,"/", date_time_sample, place_sample, "_", as.character(measured_time[i]),".csv", sep="") #the function paste is used to combine multiple character sequences. e <-try(read.csv(file_name, skip=no_skip_line, header=T), silent=FALSE) #error management if(class(e) == "try-error") next #if the file doesn't exist, skip the index; it could happen because the measurement dates may not continuous else { data_list[[i]] <- read.csv(file_name, skip=no_skip_line, header=FALSE) #read the text file and saved into the item of the list data_list[[i]] <-data_list[[i]][,-1] #Deleve the first column (A,B,...H) rownames(data_list[[i]])<-c("A", "B", "C", "D", "E", "F", "G", "H") #set names colnames(data_list[[i]])<-c("V1","V2","V3", "V4","V5","V6","V7","V8", "V9", "V10", "V11", "V12") #set names } } data_list #output (return value) of this function}エラー処理がうまくいくか確認するために以下のtime_m3には意図的に存在しないファイルの情報"34"を混ぜています。以下の3行のスクリプトを実行すれば、どのようにデータがリストに格納されるか理解できます。
#example, you can see the second element of the result list is emptytime_m3<-c(24,34,48,72,120,168)#call the function and save the result into "ecoplate_data01"ecoplate_data01<-load_ecoplate_data_osaka(folder_name="data_merged", date_time_sample="20150114-2300", place_sample ="higashihori", measured_time=time_m3, no_skip_line=7)二番目の要素[[2]]が空NULLなのがわかるでしょう。
> ecoplate_data01[[2]]NULLこれで完成ではありません。上の関数load_ecoplate_data_osakaではサンプリング日時や場所が複数ある場合には、それぞれごとに関数を呼び出してデータを読み込まなくてはいけません。比較的データ量が少ないときはそれでかまわないのですが、データ量が多い場合には、もっと自動化できることがあるはずです。そこで以下のようにcsvファイルにファイル名と対応した形のサンプリング情報の一覧を作ってみましょう。
このファイルをRスクリプトを保存しているフォルダと同じ階層の、別のフォルダ "meta_data"に保存しmeta_sample.csvと名付けたとします。まずはこのファイルを読み込みます。
#load the metadatameta_test<-read.csv("../meta_data/meta_sample.csv", header=TRUE)View(meta_test)row numberを調べれば、サンプル数も確認できます。
nrow(meta_test) #sample size in metadatameta_testにはsample_ID, data_time_sample, place_sample, measured_timeという4つの変数(4列)が格納されています。これらからファイルパスを作ることができるのですが、まずはデータのタイプを確認する必要があります。たとえば、data_time_sampleのタイプはfactorなので、
> class(meta_test$date_time_sample[2])[1] "factor"単純に数値に変換すると、factorとしての序数(1st, 2nd, 3rd等)が数値になってしまうので困ります。
> as.numeric(meta_test$date_time_sample[2])[1] 9文字列を取り出すには, as.character()で文字列に変換します。
> as.character(meta_test$date_time_sample[2])[1] "20150114-2300"以上を踏まえると以下のスクリプトでたとえば、2番目のファイルへのパス("../data_merged/20150114-2300higashihori_48.csv")を作れます。
#We can generate file pathi<-2folder_name<-"data_merged"file_ph<- paste("../" ,folder_name,"/", as.character(meta_test$date_time_sample[i]), as.character(meta_test$place_sample[i]), "_", as.character(meta_test$measured_time[i]), ".csv", sep="") file_ph以上を組み合わせて新たなファイル読み込み関数を作ってみます。
#funciton version 2 using metadata file, which stock the the metadata of results from multiple samples #data in the metadata file can be used to generate file path, instead of manually specifying as in the function above###############list of parameters########################relative_path & folder_name_data: need to specify the folder position of data #metadata: dataframe that includes the sampling date, time, and measured time information#no_skip_line: The first some lines of the raw output from microplate reader represent metadata, which should be excluded##################Output is a list; each element of the list is the data frame to have ecoplate 96 valuesload_ecoplate_data_osaka2 <- function(relative_path="../", folder_name_data="data_merged", metadata=meta_test, no_skip_line=7) { data_list <-list() #generate empty list, ecoplate data are saved as a list loop_length=nrow(metadata) #decide the loop length based on the sample size cat("The number of sample to be loaded is:", loop_length, "\n") #output basic information cat("The levels of sample places are:", levels(metadata$place_sample), "\n") #output basic information #file_name for(i in 1:loop_length) { file_name <- paste(relative_path,folder_name_data,"/", as.character(metadata$date_time_sample[i]), as.character(metadata$place_sample[i]), "_", as.character(metadata$measured_time[i]), ".csv", sep="") #the function paste is used to combine multiple character sequences. e <-try(read.csv(file_name, skip=no_skip_line), silent=FALSE) #error management if(class(e) == "try-error") next #if the file doesn't exist, skip the index; it could happen because the measurement dates may not continuous else { data_list[[i]] <- read.csv(file_name, skip=no_skip_line, header=FALSE) #read the text file and saved into the item of the list data_list[[i]] <-data_list[[i]][,-1] #Deleve the first column (A,B,...H) rownames(data_list[[i]])<-c("A", "B", "C", "D", "E", "F", "G", "H") #set names colnames(data_list[[i]])<-c("V1","V2","V3", "V4","V5","V6","V7","V8", "V9", "V10", "V11", "V12") #set names } } data_list #output (return value) of this function}関数の中の行数が長いですね。Rstudioのエディター画面では、｛｝で囲まれた領域を隠すこともできます。左端の行数が表示されている横に三角印見たいのが横向きについているのでここをマウスでクリックすると中身を隠したり再表示したいできます。
さて、このように関数を定義して、改めてメタファイルを読み込み、この関数を呼び出してosaka_dataというリストにデータを読み込んでみます。
#load metadatameta_test<-read.csv("../meta_data/meta_sample.csv", header=TRUE) #exampleosaka_data<-load_ecoplate_data_osaka2(folder_name_data="data_merged", metadata=meta_test, no_skip_line=7) 関数の実行結果はコンソールに以下のように表示されるはずです。
> osaka_data<-load_ecoplate_data_osaka2(folder_name_data="data_merged", metadata=meta_test, no_skip_line=7)The number of sample to be loaded is: 145 The levels of sample places are: higashihori honmachi 読み込まれたサンプル数が145、サンプル場所については、higashihori/honmachiの2箇所があることが分かります。
さて特定のサンプルのメタ情報とエコプレートの値を確認したければ以下のように要素をすればよいことになります。
#Check the information of metadata & list of data (osaka_data)meta_test[133,]osaka_data[[133]]コンソールには以下のように結果が表示されるでしょう。
> meta_test[133,] date_time_sample place_sample measured_time133 20150717-1100 honmachi 168> osaka_data[[133]] V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 V11 V12A 0.543 2.306 1.745 2.993 0.165 2.636 1.232 2.907 0.174 1.808 1.793 2.531B 2.674 2.927 2.010 2.844 2.625 3.194 2.285 2.771 2.036 2.615 2.237 3.201C 2.884 0.574 0.211 2.664 2.673 1.600 0.054 2.679 2.892 2.497 0.093 2.183D 2.804 1.974 2.548 2.648 2.533 2.469 2.542 2.042 2.901 2.070 2.165 2.622E 1.997 2.222 0.994 2.742 2.449 2.261 1.074 2.587 2.500 2.196 1.002 2.021F 2.466 2.060 2.656 2.161 1.803 1.781 2.655 2.279 2.468 2.635 2.564 1.907H 1.977 1.782 2.187 2.829 1.754 2.000 2.082 2.790 1.923 2.030 1.968 1.835G 2.308 1.122 0.756 2.607 1.965 1.633 2.687 1.615 2.323 1.033 2.157 1.860これでデータの読み込みは終了です。自分のサンプルの読み込みの際にはサンプル名のつけ方のルールに従って上で紹介した関数を変更すれば自動的にすべてのデータを読み込むことができるはずです。