ID轉(zhuǎn)換那點事

在進(jìn)行測序數(shù)據(jù)下游分析的時候常常需要用到不同的數(shù)據(jù)庫，而這些數(shù)據(jù)庫的分析的輸入文件經(jīng)常是各有區(qū)別胧华，因此我們常常需要將各類ID（ensamble ID、Entrez ID矩动、gene symbol等）進(jìn)行轉(zhuǎn)換。這里我們將簡單介紹使用編程進(jìn)行ID轉(zhuǎn)換的方法篮迎。

芯片ID注釋

芯片注釋主要有兩種方法：使用R包進(jìn)行注釋或在數(shù)據(jù)庫中下載芯片平臺信息然后手動注釋示姿，這兩種方法其實基本都是相同的，都是利用芯片平臺信息進(jìn)行注釋栈戳，區(qū)別只在于一個是利用R包API而另一個則是自己寫代碼。

載入數(shù)據(jù)：

> library('GEOquery')
> gse20986 <- getGEO("GSE20986", destdir=".")  # 下載處理好的表達(dá)矩陣
> exprmt <- exprs(gse20986[[1]]);head(exprmt);class(exprmt)
          GSM524662 GSM524663 GSM524664 GSM524665 GSM524666 GSM524667 GSM524668 GSM524669 GSM524670 GSM524671 GSM524672 GSM524673
1007_s_at  2.867331  3.941245  4.003509  3.230621  3.546867  2.507474  3.040758  2.456892  2.404267  3.568920  3.307409  3.225964
1053_at   10.503022  8.320770  9.523183 10.645904  8.140238 10.228986 11.019461 10.447158 10.799832  9.765222 10.069229 10.089496
117_at     2.702517  2.749144  3.226281  2.788124  4.664960  2.722543  2.745330  4.500872  2.702614  3.015139  2.768683  2.853809
121_at     3.052316  3.057630  3.056844  3.697705  3.057630  3.057630  3.057630  3.057630  3.072568  3.057630  3.057630  3.057630
1255_g_at  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998
1294_at    5.360226  4.465648  5.914568  5.558741  5.241448  5.037804  4.784192  4.978118  4.713441  4.303196  4.452989  4.766019
[1] "matrix"
# 將矩陣轉(zhuǎn)換為數(shù)據(jù)框方便后面注釋操作镊掖，并將行名稱裝換為ID
> exprmt <- as.data.frame(exprmt); exprmt$ID <- rownames(exprmt); head(exprmt)
          GSM524662 GSM524663 GSM524664 GSM524665 GSM524666 GSM524667 GSM524668 GSM524669 GSM524670 GSM524671 GSM524672 GSM524673        ID
1007_s_at  2.867331  3.941245  4.003509  3.230621  3.546867  2.507474  3.040758  2.456892  2.404267  3.568920  3.307409  3.225964 1007_s_at
1053_at   10.503022  8.320770  9.523183 10.645904  8.140238 10.228986 11.019461 10.447158 10.799832  9.765222 10.069229 10.089496   1053_at
117_at     2.702517  2.749144  3.226281  2.788124  4.664960  2.722543  2.745330  4.500872  2.702614  3.015139  2.768683  2.853809    117_at
121_at     3.052316  3.057630  3.056844  3.697705  3.057630  3.057630  3.057630  3.057630  3.072568  3.057630  3.057630  3.057630    121_at
1255_g_at  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998  2.278998 1255_g_at
1294_at    5.360226  4.465648  5.914568  5.558741  5.241448  5.037804  4.784192  4.978118  4.713441  4.303196  4.452989  4.766019   1294_at

使用annotate包注釋芯片

在數(shù)據(jù)集GSE20986頁面可以看到芯片平臺信息：

芯片平臺信息

當(dāng)然也可以直接通過代碼獲得芯片平臺信息：

> gse20986
$GSE20986_series_matrix.txt.gz
ExpressionSet (storageMode: lockedEnvironment)
assayData: 54675 features, 12 samples 
  element names: exprs 
protocolData: none
phenoData
  sampleNames: GSM524662 GSM524663 ... GSM524673 (12 total)
  varLabels: title geo_accession ... tissue:ch1 (34 total)
  varMetadata: labelDescription
featureData
  featureNames: 1007_s_at 1053_at ... AFFX-TrpnX-M_at (54675 total)
  fvarLabels: ID GB_ACC ... Gene Ontology Molecular Function (16 total)
  fvarMetadata: Column Description labelDescription
experimentData: use 'experimentData(object)'
Annotation: GPL570

知道了芯片平臺信息亩进，接下來我們就可以使用annotate包進(jìn)行注釋了。
先查看一下有沒有重復(fù)的芯片ID:

> dupID <- exprmt[duplicated(exprmt$ID),]; head(dupID); dim(dupID)
 [1] GSM524662 GSM524663 GSM524664 GSM524665 GSM524666 GSM524667 GSM524668 GSM524669 GSM524670 GSM524671 GSM524672 GSM524673 ID       
<0 行> (或0-長度的row.names)
[1]  0 13

沒有重復(fù)镐侯，直接進(jìn)行注釋

> biocLite("hgu133plus2.db")
> library(hgu133plus2.db)
> library(annotate)

# 利用getSYMBOL函數(shù)進(jìn)行注釋，這里只選擇后幾列額苟翻，方便觀察
> exprmt$symbol <- getSYMBOL(exprmt$ID, 'hgu133plus2');exprmt <- exprmt[,8:14]; head(exprmt)
          GSM524669 GSM524670 GSM524671 GSM524672 GSM524673        ID symbol
1007_s_at  2.456892  2.404267  3.568920  3.307409  3.225964 1007_s_at   <NA>
1053_at   10.447158 10.799832  9.765222 10.069229 10.089496   1053_at   RFC2
117_at     4.500872  2.702614  3.015139  2.768683  2.853809    117_at  HSPA6
121_at     3.057630  3.072568  3.057630  3.057630  3.057630    121_at   PAX8
1255_g_at  2.278998  2.278998  2.278998  2.278998  2.278998 1255_g_at GUCA1A
1294_at    4.978118  4.713441  4.303196  4.452989  4.766019   1294_at   <NA>

可以看到有一些芯片ID是沒有注釋上的崇猫，接下來對注釋情況進(jìn)行檢查需忿，主要包括三個方面：（1）未注釋上--NA诅炉；（2）未注釋上--空字符串屋厘；（3）注釋上了，但有重復(fù)议纯。

# 查看一共多少行
> dim(exprmt)
[1] 54675     7
# symbol非NA的行
> exprmt.nona <- exprmt[!is.na(exprmt$symbol),]; dim(exprmt.nona)
[1] 41941     7
> dim(exprmt[is.na(exprmt$symbol),])
[1] 12734     7
# symbol重復(fù)的行數(shù)
> symbol.dup <- exprmt.nona[duplicated(exprmt.nona$symbol),]; head(symbol.dup)
             GSM524669 GSM524670 GSM524671 GSM524672 GSM524673           ID   symbol
1552264_a_at  8.016109  7.672156  7.661139  7.835451  7.225329 1552264_a_at    MAPK1
1552272_a_at  2.279247  2.279247  2.279247  2.279247  2.279247 1552272_a_at    PRR22
1552275_s_at  4.862782  5.539661  4.672476  4.709687  3.931607 1552275_s_at      PXK
1552279_a_at  2.436048  2.731219  2.436048  2.436048  2.464805 1552279_a_at  SLC46A1
1552289_a_at  2.278998  2.278998  2.278998  2.278998  2.278998 1552289_a_at    CILP2
1552303_a_at  4.421032  4.601013  4.614977  4.628660  4.628942 1552303_a_at TMEM106A
> dim(symbol.dup)
[1] 21749     7
> exprmt.nona[exprmt.nona$symbol=='PXK',]
             GSM524669 GSM524670 GSM524671 GSM524672 GSM524673           ID symbol
1552274_at    5.290347  5.539950  5.081292  5.016670  5.032546   1552274_at    PXK
1552275_s_at  4.862782  5.539661  4.672476  4.709687  3.931607 1552275_s_at    PXK
225796_at     7.140168  7.106552  7.232615  6.650670  7.167075    225796_at    PXK
238223_at     2.281464  2.281464  2.394397  2.281464  2.281464    238223_at    PXK
# 最后看看有沒有空字符串
> exprmt.nona[exprmt.nona$symbol=='',]
[1] GSM524669 GSM524670 GSM524671 GSM524672 GSM524673 ID        symbol   
<0 行> (或0-長度的row.names)

# 最后得到非重復(fù)的溢谤，無NA的數(shù)據(jù)
> exprmt.final <- exprmt.nona[!duplicated(exprmt.nona$symbol),]; dim(exprmt.final)
[1] 20192     7
> exprmt.final[exprmt.final$symbol=='PXK',]
           GSM524669 GSM524670 GSM524671 GSM524672 GSM524673         ID symbol
1552274_at  5.290347   5.53995  5.081292   5.01667  5.032546 1552274_at    PXK

從上面可以看到，我們的芯片數(shù)據(jù)共有54675個features阀参，其中注釋到gene symbol的為41941個瞻坝，在這些注釋到的features中蛛壳，有21749個是重復(fù)的所刀。無空字符串gene symbol。最后得到非重復(fù)數(shù)據(jù)20192個（重復(fù)的則取其第一個出現(xiàn)的）赫模。

下面我們下載芯片平臺對應(yīng)信息文件手動進(jìn)行注釋：

# 下載芯片對應(yīng)信息
> gpl570 <- getGEO('GPL570', destdir=".")
> names(Meta(gpl570))
 [1] "contact_address"         "contact_city"            "contact_country"         "contact_email"           "contact_institute"       "contact_name"           
 [7] "contact_phone"           "contact_state"           "contact_web_link"        "contact_zip/postal_code" "data_row_count"          "description"            
[13] "distribution"            "geo_accession"           "last_update_date"        "manufacture_protocol"    "manufacturer"            "organism"               
[19] "relation"                "sample_id"               "series_id"               "status"                  "submission_date"         "taxid"                  
[25] "technology"              "title"                   "web_link"     
> colnames(Table(gpl570))
 [1] "ID"                               "GB_ACC"                           "SPOT_ID"                          "Species Scientific Name"         
 [5] "Annotation Date"                  "Sequence Type"                    "Sequence Source"                  "Target Description"              
 [9] "Representative Public ID"         "Gene Title"                       "Gene Symbol"                      "ENTREZ_GENE_ID"                  
[13] "RefSeq Transcript ID"             "Gene Ontology Biological Process" "Gene Ontology Cellular Component" "Gene Ontology Molecular Function"

# 選取其中的幾列用來轉(zhuǎn)換
> iddata <- Table(gpl570)[c('ID','Gene Symbol')]
> head(iddata)
         ID      Gene Symbol
1 1007_s_at DDR1 /// MIR4640
2   1053_at             RFC2
3    117_at            HSPA6
4    121_at             PAX8
5 1255_g_at           GUCA1A
6   1294_at MIR5193 /// UBA7

我們選取gpl信息中的ID蒸矛、Gene Symbol用于轉(zhuǎn)換。

> mann.anno <- merge(x=exprmt, y=iddata, by='ID',all.x=T,all.y=F); head(mann.anno); dim(mann.anno)
         ID GSM524669 GSM524670 GSM524671 GSM524672 GSM524673 symbol.RP      Gene Symbol
1 1007_s_at  2.456892  2.404267  3.568920  3.307409  3.225964      <NA> DDR1 /// MIR4640
2   1053_at 10.447158 10.799832  9.765222 10.069229 10.089496      RFC2             RFC2
3    117_at  4.500872  2.702614  3.015139  2.768683  2.853809     HSPA6            HSPA6
4    121_at  3.057630  3.072568  3.057630  3.057630  3.057630      PAX8             PAX8
5 1255_g_at  2.278998  2.278998  2.278998  2.278998  2.278998    GUCA1A           GUCA1A
6   1294_at  4.978118  4.713441  4.303196  4.452989  4.766019      <NA> MIR5193 /// UBA7

可以看到用R包注釋的（symbol.RP列）和直接下載平臺信息注釋的有些許差別斩祭，R包未能注釋到的乡话，手動注釋結(jié)果是對應(yīng)多個信息。下面我們進(jìn)一步查看一下其區(qū)別：

> library(dplyr)
> consistent <- filter(mann.anno, symbol.RP==`Gene Symbol`); head(consistent)
         ID GSM524669 GSM524670 GSM524671 GSM524672 GSM524673 symbol.RP Gene Symbol
1   1053_at 10.447158 10.799832  9.765222 10.069229 10.089496      RFC2        RFC2
2    117_at  4.500872  2.702614  3.015139  2.768683  2.853809     HSPA6       HSPA6
3    121_at  3.057630  3.072568  3.057630  3.057630  3.057630      PAX8        PAX8
4 1255_g_at  2.278998  2.278998  2.278998  2.278998  2.278998    GUCA1A      GUCA1A
5   1316_at  5.246463  5.547843  5.269052  4.794183  4.620719      THRA        THRA
6   1320_at  4.578203  4.446357  4.492585  4.502381  4.446357    PTPN21      PTPN21
> dim(mann.anno); dim(consistent)
[1] 54675     8
[1] 39145     8
# 可以看到一共有39145行注釋結(jié)果都是相同的，那么這些中有多少是都未注釋到的呢（即NA）
> consistent.na <- filter(consistent, is.na(symbol.RP));head(consistent.na)
[1] ID          GSM524669   GSM524670   GSM524671   GSM524672   GSM524673   symbol.RP   Gene Symbol
<0 行> (或0-長度的row.names)
# 結(jié)果是這些相同的中沒有NA
# 再來看看注釋到重復(fù)gene的情況吧
> dup <- consistent[duplicated(consistent$symbol.RP),]
> head(dup);dim(dup)
             ID GSM524669 GSM524670 GSM524671 GSM524672 GSM524673 symbol.RP Gene Symbol
16 1552264_a_at  8.016109  7.672156  7.661139  7.835451  7.225329     MAPK1       MAPK1
20 1552272_a_at  2.279247  2.279247  2.279247  2.279247  2.279247     PRR22       PRR22
22 1552275_s_at  4.862782  5.539661  4.672476  4.709687  3.931607       PXK         PXK
26 1552279_a_at  2.436048  2.731219  2.436048  2.436048  2.464805   SLC46A1     SLC46A1
32 1552289_a_at  2.278998  2.278998  2.278998  2.278998  2.278998     CILP2       CILP2
40 1552303_a_at  4.421032  4.601013  4.614977  4.628660  4.628942  TMEM106A    TMEM106A
[1] 20267     8
# 因此最后非重復(fù)的
> dim(consistent[!duplicated(consistent$symbol.RP),])
[1] 18878     8

對兩種注釋結(jié)果相同的 gene symbol 進(jìn)行分析屋群，一共得到了39145個注釋上的gene坏挠，這其中有20267個是重復(fù)的，因此最后得到了18878個非重復(fù)注釋的gene降狠。
我們最后再對手動注釋的結(jié)果進(jìn)行分析一下。

# 看看為注釋到的情況（NA）
> mann.anno[is.na(mann.anno$`Gene Symbol`),]
[1] ID          GSM524669   GSM524670   GSM524671   GSM524672   GSM524673   symbol.RP   Gene Symbol
<0 行> (或0-長度的row.names)
# 沒有注釋不上的否纬，看來R包注釋不上的這里都注釋到多個gene symbol了（如DDR1 /// MIR4640）
> library(stringr)
> multi.anno <- mann.anno[str_detect(mann.anno$`Gene Symbol`, '///'),]; head(multi.anno); dim(multi.anno)
              ID GSM524669 GSM524670 GSM524671 GSM524672 GSM524673 symbol.RP                Gene Symbol
1      1007_s_at  2.456892  2.404267  3.568920  3.307409  3.225964      <NA>           DDR1 /// MIR4640
6        1294_at  4.978118  4.713441  4.303196  4.452989  4.766019      <NA>           MIR5193 /// UBA7
16    1552258_at  2.437344  2.458223  2.437344  2.437344  2.429044     CYTOR LINC00152 /// LOC101930489
32  1552283_s_at  2.278998  2.278998  2.278998  2.281762  2.278998      <NA>       ZDHHC11 /// ZDHHC11B
97    1552388_at  2.550713  2.389557  2.389557  2.389557  2.389557      <NA>        FLJ30901 /// SCUBE1
108 1552401_a_at  2.278998  2.278998  2.278998  2.294001  2.278998      <NA>        BACH1 /// GRIK1-AS2
[1] 2796    8
# 這里2796個features注釋到多gene symbol蛋褥，但意外的是，有的多gene symbol位點R包也注釋到了
# 看看dup情況怎么樣
> mann.anno.rmdup <- mann.anno[!duplicated(mann.anno$`Gene Symbol`),]; head(mann.anno.rmdup); dim(mann.anno.rmdup)
         ID GSM524669 GSM524670 GSM524671 GSM524672 GSM524673 symbol.RP      Gene Symbol
1 1007_s_at  2.456892  2.404267  3.568920  3.307409  3.225964      <NA> DDR1 /// MIR4640
2   1053_at 10.447158 10.799832  9.765222 10.069229 10.089496      RFC2             RFC2
3    117_at  4.500872  2.702614  3.015139  2.768683  2.853809     HSPA6            HSPA6
4    121_at  3.057630  3.072568  3.057630  3.057630  3.057630      PAX8             PAX8
5 1255_g_at  2.278998  2.278998  2.278998  2.278998  2.278998    GUCA1A           GUCA1A
6   1294_at  4.978118  4.713441  4.303196  4.452989  4.766019      <NA> MIR5193 /// UBA7
[1] 23521     8

從上面可以看到所有芯片ID都注釋到了谬俄，沒有NA，其中沒有重復(fù)的有23521個（其中2796個是多gene symbol溃论，因此單個gene symbol 20725個）痘昌，較R包注釋到的（20276）要多一些钥勋。

ID轉(zhuǎn)換

接下來我們將上面注釋到的Gene symbol轉(zhuǎn)換為Entrez ID辆苔，使用的是HGNC提供的ID對應(yīng)文件。

> hgnc <- read.csv('D:/TCGA/hgnc_complete_set.txt',header = T,stringsAsFactors = F,sep='\t');colnames(hgnc)
 [1] "hgnc_id"                  "symbol"                   "name"                     "locus_group"              "locus_type"               "status"                  
 [7] "location"                 "location_sortable"        "alias_symbol"             "alias_name"               "prev_symbol"              "prev_name"               
[13] "gene_family"              "gene_family_id"           "date_approved_reserved"   "date_symbol_changed"      "date_name_changed"        "date_modified"           
[19] "entrez_id"                "ensembl_gene_id"          "vega_id"                  "ucsc_id"                  "ena"                      "refseq_accession"        
[25] "ccds_id"                  "uniprot_ids"              "pubmed_id"                "mgd_id"                   "rgd_id"                   "lsdb"                    
[31] "cosmic"                   "omim_id"                  "mirbase"                  "homeodb"                  "snornabase"               "bioparadigms_slc"        
[37] "orphanet"                 "pseudogene.org"           "horde_id"                 "merops"                   "imgt"                     "iuphar"                  
[43] "kznf_gene_catalog"        "mamit.trnadb"             "cd"                       "lncrnadb"                 "enzyme_id"                "intermediate_filament_db"
[49] "rna_central_ids"     

> destid <- hgnc[c('symbol','entrez_id','ensembl_gene_id')];head(destid)
    symbol entrez_id ensembl_gene_id
1     A1BG         1 ENSG00000121410
2 A1BG-AS1    503538 ENSG00000268895
3     A1CF     29974 ENSG00000148584
4      A2M         2 ENSG00000175899
5  A2M-AS1    144571 ENSG00000245105
6    A2ML1    144568 ENSG00000166535

# ID 轉(zhuǎn)換
> trans <- merge(x=exprmt[!duplicated(exprmt$symbol.RP),], y=destid, by.x='symbol.RP',by.y='symbol');cat('exprmt.rmdup\n');dim(exprmt[!duplicated(exprmt$symbol.RP),]);cat('destid\n');dim(destid);cat('trans\n');dim(trans)
exprmt.rmdup
[1] 20193     7
destid
[1] 42508     3
trans
[1] 19062     9

> trans.rmdup <- trans[!duplicated(trans$symbol.RP),]; dim(trans.rmdup)
[1] 19061     9

> trans.rmdup <- trans[!duplicated(trans$entrez_id),]; dim(trans.rmdup)
[1] 19060     9

可以看到，使用HGNC提供的ID轉(zhuǎn)換文件骑冗，我們共得到了19062個無重復(fù)ID，比上面注釋到的少了快2000個贼涩，對于這些差異基因有機(jī)會在探究一下。

最后編輯于：2018.05.12 15:33:00

?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者

人面猴
序言：七十年代末谤绳，一起剝皮案震驚了整個濱河市，隨后出現(xiàn)的幾起案子缩筛，更是在濱河造成了極大的恐慌，老刑警劉巖歪脏，帶你破解...
沈念sama閱讀 219,490評論 6贊 508
死咒
序言：濱河連續(xù)發(fā)生了三起死亡事件，死亡現(xiàn)場離奇詭異，居然都是意外死亡啄寡，警方通過查閱死者的電腦和手機(jī)，發(fā)現(xiàn)死者居然都...
沈念sama閱讀 93,581評論 3贊 395
救了他兩次的神仙讓他今天三更去死
文/潘曉璐我一進(jìn)店門挺物，熙熙樓的掌柜王于貴愁眉苦臉地迎上來，“玉大人砚著，你說我怎么就攤上這事』拢” “怎么了？”我有些...
開封第一講書人閱讀 165,830評論 0贊 356
道士緝兇錄：失蹤的賣姜人
文/不壞的土叔我叫張陵舌镶，是天一觀的道長。經(jīng)常有香客問我餐胀，道長，這世上最難降的妖魔是什么否灾？我笑而不...
開封第一講書人閱讀 58,957評論 1贊 295
?港島之戀（遺憾婚禮）
正文為了忘掉前任鸣奔，我火速辦了婚禮墨技，結(jié)果婚禮上溃蔫，老公的妹妹穿的比我還像新娘。我一直安慰自己伟叛，他們只是感情好，可當(dāng)我...
茶點故事閱讀 67,974評論 6贊 393
惡毒庶女頂嫁案：這布局不是一般人想出來的
文/花漫我一把揭開白布紊遵。她就那樣靜靜地躺著账千，像睡著了一般暗膜。火紅的嫁衣襯著肌膚如雪。梳的紋絲不亂的頭發(fā)上学搜，一...
開封第一講書人閱讀 51,754評論 1贊 307
城市分裂傳說
那天，我揣著相機(jī)與錄音聚磺，去河邊找鬼炬丸。笑死瘫寝，一個胖子當(dāng)著我的面吹牛稠炬，可吹牛的內(nèi)容都是我干的。我是一名探鬼主播首启，決...
沈念sama閱讀 40,464評論 3贊 420
雙鴛鴦連環(huán)套：你想象不到人心有多黑
文/蒼蘭香墨我猛地睜開眼，長吁一口氣：“原來是場噩夢啊……” “哼栽惶！你這毒婦竟也來了？” 一聲冷哼從身側(cè)響起外厂，我...
開封第一講書人閱讀 39,357評論 0贊 276
萬榮殺人案實錄
序言：老撾萬榮一對情侶失蹤代承，失蹤者是張志新（化名）和其女友劉穎，沒想到半個月后论悴，有當(dāng)?shù)厝嗽跇淞掷锇l(fā)現(xiàn)了一具尸體，經(jīng)...
沈念sama閱讀 45,847評論 1贊 317
?護(hù)林員之死
正文獨居荒郊野嶺守林人離奇死亡膀估，尸身上長有42處帶血的膿包…… 初始之章·張勛以下內(nèi)容為張勛視角年9月15日...
茶點故事閱讀 37,995評論 3贊 338
?白月光啟示錄
正文我和宋清朗相戀三年察纯，在試婚紗的時候發(fā)現(xiàn)自己被綠了帕棉。大學(xué)時的朋友給我發(fā)了我未婚夫和他白月光在一起吃飯的照片。...
茶點故事閱讀 40,137評論 1贊 351
活死人
序言：一個原本活蹦亂跳的男人離奇死亡慰枕，死狀恐怖，靈堂內(nèi)的尸體忽然破棺而出具帮，到底是詐尸還是另有隱情低斋，我是刑警寧澤蜂厅，帶...
沈念sama閱讀 35,819評論 5贊 346
?日本核電站爆炸內(nèi)幕
正文年R本政府宣布膊畴，位于F島的核電站，受9級特大地震影響巴比，放射性物質(zhì)發(fā)生泄漏礁遵。R本人自食惡果不足惜，卻給世界環(huán)境...
茶點故事閱讀 41,482評論 3贊 331
男人毒藥：我在死后第九天來索命
文/蒙蒙一政勃、第九天我趴在偏房一處隱蔽的房頂上張望。院中可真熱鬧奸远，春花似錦、人聲如沸懒叛。這莊子的主人今日做“春日...
開封第一講書人閱讀 32,023評論 0贊 22
一樁弒父案薛窥，背后竟有這般陰謀
文/蒼蘭香墨我抬頭看了看天上的太陽眼姐。三九已至，卻和暖如春众旗，著一層夾襖步出監(jiān)牢的瞬間，已是汗流浹背贡歧。一陣腳步聲響...
開封第一講書人閱讀 33,149評論 1贊 272
情欲美人皮
我被黑心中介騙來泰國打工，沒想到剛下飛機(jī)就差點兒被人妖公主榨干…… 1. 我叫王不留持际，地道東北人。一個月前我還...
沈念sama閱讀 48,409評論 3贊 373
代替公主和親
正文我出身青樓蜘欲，卻偏偏與公主長得像，于是被迫代替她去往敵國和親姥份。傳聞我的和親對象是個殘疾皇子，可洞房花燭夜當(dāng)晚...
茶點故事閱讀 45,086評論 2贊 355