ggplot2でdrcモデルをプロットする。 seq（）との問題

Question

私のモデルは、Rベースのグラフィックスではありますが、ggplot2でプロットしたときに漸近線に向かって進まない。 ggplot2では、X軸上の特定の点（画像は以下を含む）で停止し、私は、それが私が変換するpredict()を使用しています。このポスト

の底部にポストされたデータをこれがseq()に関連している90％確信drm（用量応答パッケージ）ロジットモデルからのデータ。基本グラフィックスでは、ggplot2ではS字状のカーブがすばらしく見えます。
（私はあなたがそれを必要としないと思うが、念のために、データは、この記事の一番下に含まれている）

library(drc) 
library(ggplot2) 

ロジットモデルにデータをフィット

mod1 <- drm(probability ~ (dose), weights = total, data = mydata1,  type ="binomial", fct = LL.2()) 

plot(mod1,broken=FALSE,type="all",add=FALSE, col= "purple", xlim=c(0, 10000)) 

Image of Base graphic 2-parameter logit

著者のデモのコードを使用してモデルのデータを抽出しています（下記リンク先）：

newdata1 <-expand.grid(dose=exp(seq(log(0.5),log(100),length=200))) 

pm1<- predict(mod1, newdata=newdata1,interval="confidence") 

newdata1$p1 <-pm1[,1] 

newdata1$pmin1 <-pm1[,2] 

newdata1$pmax1 <- pm1[,3] 

最後にggplot2グラフィック。

p1 <- ggplot(mydata1, aes(x=dose01,y=probability))+ 
    geom_point()+ 
    geom_ribbon(data=newdata1, aes(x=dose,y=p1, 
ymin=pmin1,ymax=pmax1),alpha=0.2,color="blue",fill="pink") + 
    geom_step(data=newdata1, aes(x=dose,y=p1))+ 
    coord_trans(x="log") + #creates logline for x axis 
    xlab("dose")+ylab("response") 

Images 1&2 and 3&4 seqによって私のプロットの違いを示します。

seq()に以下を使用すると、グラフィックがデータよりも先に引き出されます。 (seq(log(0.5),**log(10000)**,length=200)))（画像1 & 2）

私はそれを研究にもかかわらず、seq()を理解していません。誰かが私のプロットで起こっていることを理解するのに役立ちますか？

seqの最初の用語は下限を定義しているようです。 3番目の用語は何を定義していますか？あなたは画像でこれを見ることができます3 & 4 - グラフィックはまともではありません。私はこの問題をやや難解にしてきましたが、それでもなおinfinに向かって進まないのです。私は8つのロジットモデルを共同作図するので、これはわずかな問題です。

- トラブルggplot2でDRC/DRMモデルをプロットした人のため （私はこれらのタイトルを持つので、裸つ以上のリンクを投稿することはできません）

、以下の記事は非常に有用だった：検索 プロット用量 - 応答曲線-と-ggplot2-と-DRC
と、このタイトルのために：プロット用量 - 応答曲線-と-ggplot2-と-DRC

私が続いてきましたDRCの著者は、彼の記事のサポート情報に記載されています - このコードの一部は上記で使用されています。論文タイトル：R、Christopher Ritzを用いた線量 - 応答解析PlosOne。

このデータは理想的な形式よりも少ない、私は

> dput(mydata1) 

structure(list(dose = c(25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 75L, 75L, 75L, 75L, 75L, 75L, 
75L, 75L, 75L, 75L, 75L, 75L, 75L, 75L, 75L, 100L, 100L, 100L, 
100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
100L, 150L, 150L, 150L, 150L, 150L, 150L, 150L, 150L, 150L, 150L, 
150L, 150L, 150L, 150L, 150L, 200L, 200L, 200L, 200L, 200L, 200L, 
200L, 200L, 200L, 200L, 200L, 200L, 200L, 200L, 200L), total = c(25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L), affected = c(2L, 
0L, 0L, 0L, 0L, 1L, 1L, 2L, 2L, 4L, 1L, 1L, 4L, 0L, 10L, 0L, 
1L, 0L, 1L, 0L, 3L, 0L, 4L, 2L, 0L, 2L, 0L, 1L, 2L, 3L, 2L, 0L, 
2L, 0L, 0L, 4L, 0L, 1L, 2L, 3L, 0L, 21L, 1L, 3L, 1L, 2L, 7L, 
0L, 0L, 0L, 0L, 8L, 7L, 3L, 7L, 2L, 2L, 10L, 3L, 4L, 0L, 7L, 
0L, 3L, 3L, 20L, 25L, 22L, 23L, 22L, 18L, 14L, 20L, 20L, 21L), 
    probability = c(0.08, 0, 0, 0, 0, 0.04, 0.04, 0.08, 0.08, 
    0.16, 0.04, 0.04, 0.16, 0, 0.4, 0, 0.04, 0, 0.04, 0, 0.12, 
    0, 0.16, 0.08, 0, 0.08, 0, 0.04, 0.08, 0.12, 0.08, 0, 0.08, 
    0, 0, 0.16, 0, 0.04, 0.08, 0.12, 0, 0.84, 0.04, 0.12, 0.04, 
    0.08, 0.28, 0, 0, 0, 0, 0.32, 0.28, 0.12, 0.28, 0.08, 0.08, 
    0.4, 0.12, 0.16, 0, 0.28, 0, 0.12, 0.12, 0.8, 1, 0.88, 0.92, 
    0.88, 0.72, 0.56, 0.8, 0.8, 0.84)), .Names = c("dose", "total", 
"affected", "probability"), row.names = c(NA, -75L), class = "data.frame") 

Answer 1

これは長いコメントです。私はあなたが間違っての値を間違えたと思うpredict()またはaes(x)の値を与える。

log10000 <- exp(seq(log(0.5), log(10000), length=200)) 
log1000 <- exp(seq(log(0.5), log(1000), length=200)) 

log10000df <- as.data.frame(cbind(dose = log10000, predict(mod1, data.frame(dose = log10000), interval="confidence"))) 
log1000df <- as.data.frame(cbind(dose = log1000, predict(mod1, data.frame(dose = log1000), interval="confidence"))) 

## a common part 
p0 <- ggplot(mydata1, aes(x = dose, y = probability)) + 
    geom_point() + coord_trans(x="log") + 
    xlab("dose") + ylab("response") + xlim(0.5, 10001) 

p10000 <- p0 + geom_line(data = log10000df, aes(x = dose, y = Prediction)) + 
    geom_ribbon(data = log10000df, aes(x = dose, y = Prediction, ymin = Lower, ymax = Upper), 
       alpha = 0.2, color = "blue", fill = "pink") 

p1000 <- p0 + geom_line(data = log1000df, aes(x = dose, y = Prediction)) + 
    geom_ribbon(data = log1000df, aes(x = dose, y = Prediction, ymin = Lower, ymax = Upper), 
       alpha = 0.2, color = "blue", fill = "pink") 

Answer 2

をknow-せている場合は、seqの用語が何をするかの説明については?seqを参照してください。 seq(log(.5), log(1000), length=200)は、log（0.5）からlog（1000）になる200の数値を均等に配置します。 3番目の引数を指定しない場合（またはby=XYZを指定した場合）、数値間の間隔です。

したがって、seq(log(0.5), log(1000), length=200)を実行すると、log（0.5）からlog（1000）までの200ポイントで適合を計算します。

「無限に行く」とは、リンクされた画像のように、エッジの前で停止するのではなく、プロットの端からラインを外すことを意味します。デフォルトでは、ggplotはプロットされたすべてのものがプロットに収まるようにしようとするので、データの範囲を少し超えて軸を拡張します。

制限する場合は、+ xlim(c(lower, upper))を使用してください。

（私はあなたの例を再現することはできませんのでdrcをインストールするトラブルを抱えています。ここでは、おもちゃ一つだ）上記

x = seq(0.5, 100, length=200) 
df <- data.frame(x=x, y=x^2) 
ggplot(df, aes(x=x, y=y)) + geom_line() + coord_trans(x="log") 

、ラインは（予想通り）100とに出て行きます限界は少し超えています。私は、プロットの端に触れるためにラインを望んでいた場合、私は（例えば）ちょうど100でのxの上限をクリップすることができます - coord_transにlimx引数を使用します。

ggplot(df, aes(x=x, y=y)) + geom_line() + coord_trans(x="log", limx=c(0.5, 100)) 

ですから、プロットするときあなたのモデルは、あなたのx軸の境界線がどのようなものかを決定し、それらの値に沿ってすべてのモデルを予測するようにしてください。そして、xの限界をそれらの境界に制限すると、線は「無限に進む」ように見えます。