2016-04-05 25 views
2

私のロジスティック回帰モデル、すなわち「S」形状のロジスティック曲線をプロットするためのコードを記述します。 2つの独立した共変量を持っているので、どうすればいいですか?データセットとモデルのコードを添付しています。前もって感謝します。ロジスティック回帰モデルのロジスティック曲線を作成する

239 0.72 1 
324.6 0.83 1 
331.8 0.95 1 
334.3 0.83 1 
259.7 0.89 1 
212.3 0.88 1 
204.7 0.65 1 
253.86 0.75 1 
258.94 0.85 1 
329.66 0.95 0 
469.68 1.46 0 
459.74 1.11 0 
293.2 0.64 0 
297.88 0.98 0 
267.9 0.82 0 
374.1 1.29 0 
333.62 0.74 0 


dat <- read.table("data.txt") 
colnames(dat)<-c("press","v","gender") 

# logostic regression 
dat$gender <- factor(dat$gender) 
mylogit<- glm(gender~press+v,data=dat,family="binomial") 
summary(mylogit) 

######## the code below are irrelevant to making plot, ignore if you want 

mylogit$fitted.values 

newdat <- data.frame(t(c(300,0.1))) 
colnames(newdat)<-c("press","v") 
    # this is your new dataset, we name it as "newdat" 
pred <- predict(mylogit,newdata = newdat,type="response") 
pred # the probability of being in class 1 will stored in this object 

pred <- predict(mylogit,newdata = dat,type="response") 
pred # the probability of being in class 1 will stored in this object 
# accuracy 
dat$pred <- 0 
factor(dat$pred) 
dat$pred[which(pred>0.5)] <- 1 

table(dat$gender,dat$pred) 
+0

はい - ロジスティック曲線をプロットしたい –

+0

これはこのようなものでなければなりませんが、エラーが表示されています。私はそれをどう扱うかわかりません。 –

+0

plot(dat $ press、dat $ gender) カーブ(予測(mylogit、newdata = newdat、type = "resp")、add = TRUE) –

答えて

10

あなたは2つの連続した非カテゴリ変数を持っているので、ロジスティック曲線は3D曲線になります。私はあなたに提示のための2つの方法を提供します。

  • persp実際の3D滑らかな曲線を生成する機能。
  • vを値の数に固定し、2Dロジスティック曲線(これを "S"字型曲線と呼びます)を生成します。

あなたが最初の2Dグリッドを生成する必要が3D曲線

press_grid <- seq(200, 480, by = 5) 
v_grid <- seq(0.6, 1.5, by = 0.1) 
newdat <- data.frame(press = rep(press_grid, times = length(v_grid)), v = rep(v_grid, each = length(press_grid))) 
pred <- predict.glm(mylogit, newdata = newdat, type="response") 
z <- matrix(pred, length(press_grid)) 
persp(press_grid, v_grid, z, xlab = "pressure", ylab = "velocity", zlab = "predicted probability", main = "logistic curve (3D)", theta = 30, phi = 20) 

newdatにはこのグリッドがあり、このグリッドを表示するにはplot(newdat)を実行できます。次に、predict.glm(..., type = "response")を呼び出すことによって、このグリッド上で予測が行われます。結果はpredです。それをプロットするには、マトリックスzにキャストし、次にperspを呼び出して3Dプロットを作成します。 xlab,ylabおよびzlabは、3軸のラベルである。パラメータthetaphiは、視野角を調整するために使用されます。上記で

pressv用限界グリッドは元のデータの範囲に基づいている:range(dat$press)range(dat$v)。この範囲を超えた予測はあまりにも遠すぎません。しかし、この範囲内でさえ、あなたは17の観察しかありません。だからあなたはまだプロットに懐疑的である必要があります。ここ

は曲線である:

3D

2D曲線

このおもちゃの機能は、いくつかのレベルに固定vで、2次元曲線を作製するために有用である:

curve_2D_fix_v <- function(model, v = 1, press_grid = seq(200, 480, by = 5), add = FALSE, col = "black") { 
    newdat <- data.frame(press = press_grid, v = v) 
    pred <- predict.glm(model, newdat, type = "response") 
    if (add) lines(press_grid, pred, col = col) else { 
    plot(press_grid, pred, xlab = "pressure", ylab = "predicted probability", type = "l", col = col, main = "logistic curve (2D)") 
    abline(h = c(0, 0.5, 1), lty = 2, col = col) 
    } 
    } 

add = FALSEの場合、新しいプロットウィンドウが開きます。 TRUEですが、それは前のウィンドウにプロットします(しかし、そのようなウィンドウがあることを確認するのがあなたの義務です!)2Dプロットは、0、0.5、1で水平線を追加できるので、より多くの情報を与えます。

行こてみましょう。ここでは

curve_2D_fix_v(mylogit, v = 0.4, add = FALSE, col = "black") 
curve_2D_fix_v(mylogit, v = 0.6, add = TRUE, col = "red") 
curve_2D_fix_v(mylogit, v = 0.8, add = TRUE, col = "green") 
curve_2D_fix_v(mylogit, v = 1, add = TRUE, col = "blue") 
curve_2D_fix_v(mylogit, v = 1.2, add = TRUE, col = "cyan") 
curve_2D_fix_v(mylogit, v = 0.4, add = TRUE, col = "yellow") 

が曲線である:両方のプロットで

2D


ディスカッション

、我々は間の関係を見ますgender(予測確率)およびv(速度)はあまり強くありません。 2Dプロットでは、vのほぼすべての値が同じ曲線を生成します。一方、press(圧力)は強い効果である。戻るあなたのモデルに

:あなたはvは全く重要ではないことがわかります

> summary(mylogit) 
Coefficients: 
      Estimate Std. Error z value Pr(>|z|) 
(Intercept) 8.08326 4.45463 1.815 0.0696 . 
press  -0.02575 0.01618 -1.591 0.1115 
v   -0.15385 4.83824 -0.032 0.9746 
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

!厳密に言えば、pressも0.1レベルで有意ではない。 これは非常に弱いモデルです。変数vをドロップして、pressを唯一の変数として使用してモデルを再度実行することをお勧めします。

+0

有益な情報をありがとうございました。ロジスティック回帰は数値変数** press **でしか動かないとは思っていませんでしたが、私はあなたの提案に従います。 –

関連する問題