・お題：何らかの反応の阻害に関して、96穴ベースのアッセイで結果を得た。50%阻害濃度などを求めてみたい。
・なお、当方は専門家でもなんでもないので、解析が正しいかは分からない。
・適当にデータを作った。今回はEXCELフォーマットを読み込んで試験物質の名称、濃度、アッセイスコアを読み込み、各試験物質の濃度とアッセイスコアの関係からDose Response Curveを作成する。
・試験物質情報
・濃度情報
・アッセイスコア
・ファイルをカレントディレクトリに放り込み、データを読み込む。
> library(tidyverse)
> library(openxlsx)
> xlsx_file = "DrcSample.xlsx"
> Compound =
+ read.xlsx(xlsxFile = xlsx_file,
+ sheet = "Compound",
+ rowNames = TRUE) %>%
+ pivot_longer(data = .,
+ cols = colnames(.),
+ names_to = "colm",
+ values_to = "value")
>
> Concentration =
+ read.xlsx(xlsxFile = xlsx_file,
+ sheet = "Concentration",
+ rowNames = TRUE) %>%
+ pivot_longer(data = .,
+ cols = colnames(.),
+ names_to = "colm",
+ values_to = "value")
>
> Score =
+ read.xlsx(xlsxFile = xlsx_file,
+ sheet = "Score",
+ rowNames = TRUE) %>%
+ pivot_longer(data = .,
+ cols = colnames(.),
+ names_to = "colm",
+ values_to = "value")
>
> df =
+ data.frame(Compound = Compound$value,
+ Concentration = Concentration$value,
+ Score = Score$value)
> df %>% head(.)
・理想的には0-100%のシグモイダルカーブになると分かっているとする。この時、0%と100%の値を平均値で算出しておき、後で外挿する。
> df_mean = df %>%
+ group_by(., Compound) %>%
+ summarise_at(., vars(Score), list(mean = mean))
> df_mean %>% head()
> ctrl_0 = df_mean %>%
+ filter(., Compound == "0_ctrl") %>%
+ .$mean
> ctrl_100 = df_mean %>%
+ filter(., Compound == "100_ctrl") %>%
+ .$mean
・次に、試験物質の一覧を得ておく。
> Compound_List = df$Compound %>% unique(.)
> Compound_List = Compound_List[- which(Compound_List %in% c("100_ctrl", "0_ctrl"))]
> Compound_List
[1] "A001" "A003" "A002" "A004"
・ED50を一気に算出して表にまとめるなら以下のような感じ。drcパッケージを使う。
> library(drc)
> df_ED50 = NULL
> for (variable in Compound_List) {
+ df_target = df %>% dplyr::filter(., Compound == variable)
+ LL4 <- drm(Score ~ Concentration, #Log 4para LogisticsでFitting
+ data = df_target,
+ fct = LL.4(fixed = c(NA, ctrl_0, ctrl_100, NA), #ここで外挿するか否かはケースバイケース。今回はやり方をメモしておきたいのでやる版
+ names = c("b", "c", "d", "e"))) #c:min, d:max, e:ed50
+ df_ED50 = rbind(df_ED50, ED(object = LL4,
+ respLev = c(50),
+ interval = "delta",
+ level = 0.95))
+ }
> rownames(df_ED50) = Compound_List
> df_ED50
・適当にプロットしてみる。
> for (variable in Compound_List) {
+ df_target = df %>% dplyr::filter(., Compound == variable)
+ LL4 <- drm(Score ~ Concentration,
+ data = df_target,
+ fct = LL.4(fixed = c(NA, ctrl_0, ctrl_100, NA),
+ names = c("b", "c", "d", "e")))
+ plot(LL4,
+ main = variable,
+ type = "all",
+ ylim = c(ctrl_0 - 10, ctrl_100 + 10))
+ plot(LL4,
+ type="confidence",
+ col = "lightblue",
+ add=TRUE)
+ }
・なんだかそれっぽい結果になったっぽい。
おわり
