2016-04-19 6 views
1

これは非常にばかげた質問かもしれないと思うが、ここに行く! (CrossValidatedにもっと適している場合は申し訳ありませんが、これはプログラミング上の問題であるか、より統計的に近いものであるかはわかりません...ステップ関数を90度の遷移にする

私はcumSeg階段関数、別名、小定数定数関数)を使用して、下のコード/図のように不連続な(x軸の)データにフィットさせてください。

それはすべて問題ありません。階段関数(赤い線)に垂直な遷移を持たせることができます(つまり、関数の両方の肩が90度になるようにすることができます)。これを行うには、x軸の値が現在の2ジャンプポイント。これは可能ですか?

もしそうなら、私は別の質問に繋がります。これらの90度の遷移と垂直降下がある場合、チャート上のその線のずれ?

# Plotting step-functions on to GC-operon data. 

require(ggplot2) 
library("ggplot2") 
require(reshape2) 
library("reshape2") 
require(scales) 
library(RColorBrewer) 
library(cumSeg) 

df <- structure(list(PVC1 = 0.4019026, PVC2 = 0.4479259, PVC3 = 0.4494118, PVC4 = 0.4729437, 
       PVC5 = 0.4800556, PVC6 = 0.449229, PVC7 = 0.4905295, PVC8 = 0.4457566, 
       PVC9 = 0.4271259, PVC10 = 0.4850341, PVC11 = 0.4369965, PVC12 = 0.4064052, 
       PVC13 = 0.3743776, PVC14 = 0.3603853, PVC15 = 0.3965469, PVC16 = 0.365461), 
       .Names = c("PVC1","PVC2","PVC3","PVC4","PVC5","PVC6","PVC7","PVC8", 
          "PVC9","PVC10","PVC11","PVC12","PVC13","PVC14","PVC15","PVC16"), 
       class = "data.frame", 
       row.names = c(NA, -1L) 
      ) 

melted_df <- melt(df, variable.name = "Locus", value.name = "GC") 
st_dev <- c(0.023031363, 0.024919217, 0.017371129, 
     0.019008759, 0.026650605, 0.026904926, 
     0.024227542, 0.017767553, 0.026152478, 
     0.039770898, 0.023929714, 0.028845442, 
     0.015572219, 0.024967336, 0.014955416, 0.024569096) 

operon_gc <- 0.408891366 
opgc_stdev <- 0.015712091 
genome_gc <- 0.425031611 
gengc_stdev <- 0.007587437 

stepfunc <- jumpoints(y=melted_df$GC*100, k=1, output="1") 

gc_chart <- ggplot(melted_df, aes(Locus, GC*100, fill=Locus,)) + 
       geom_bar(width=0.6, stat = "identity") 
gc_chart <- gc_chart + ylab("GC Content (%)") 
gc_chart <- gc_chart + theme(axis.text.x = element_text(angle=45, hjust=1)) 
gc_chart <- gc_chart + geom_abline(intercept=operon_gc*100, 
           slope=0, 
           colour="gray", 
           linetype=3, 
           show.legend =TRUE) 
gc_chart <- gc_chart + geom_text(aes(15.7, 41.7, label="Operon GC"), 
          size=5, 
          color="gray") 
gc_chart <- gc_chart + geom_abline(intercept=genome_gc*100, 
           slope=0, 
           colour="black", 
           linetype=3, 
           show.legend = TRUE) 
gc_chart <- gc_chart + geom_text(aes(15.7, 43.3, label="Genome GC"), 
          size=5, 
          color="black") 
gc_chart <- gc_chart + coord_cartesian(ylim=c(30,55)) 
gc_chart <- gc_chart + geom_errorbar(width=.2, size=0.4, color="azure4", aes(Locus, 
       ymin = (GC - cbind(melted_df, st_dev)$st_dev)*100, 
       ymax = (GC + cbind(melted_df, st_dev)$st_dev)*100)) 
gc_chart <- gc_chart + geom_line(linetype=2,aes(x=as.numeric(Locus), y=stepfunc$fitted.values, colour="red", group=1)) 

gc_chart 

enter image description here

EDIT:@Gregor、geom_stepは文字通り私のコードでstepためlineを代入(ありがとう、この生成所望の効果を達成した。しかし enter image description here

を、それグラフ、ブレークポイントはPVC12になりますが、関数からブレークポイント値を抽出するとき...

> stepfunc$psi 
V 
11 

この場合、グラフが誤解を招く恐れがあります。おそらく、11と12の間にブレークがあることを示す以前のバージョンでは、私はもっと良いでしょう。

+1

'?geom_step'を参照してください。 – Gregor

+0

ブレイクポイントの考えは+1をシフトしますか? –

+0

答えを見てください。また、 'library()'は同じパッケージを必要とする 'require()'が冗長です。最も一般的なアドバイスは常に 'library()'を使うことです。 – Gregor

答えて

1

ggplotは、データをプロットするのに最適です。 xの値は1:12で、yの値は45と37の値です。 yの値(およびy値の変更)が整数値のxで定義されていない場合は、値を変更してください。

step_df = data.frame(x = c(1, 11.5, 16), y = c(45.24568, 37.5, 37.5)) 

gc_chart + geom_step(data = stepdf, linetype=2, 
        aes(x = x, y = y, colour="red", group=1), 
        inherit.aes = F) 

私はプログラムで適切なstep_dfを定義しておきます。