作业要求：

根据给定的两个数据文件genefile.txt（基因表达谱）和geneid（基因ID号），利用熵值评价各个基因的表达稳定性，并找出前10个表达最不稳定的基因（可能差异表达）。

基本上就是对熵的理解，和R中entropy包的使用
↓该加载加载，该读数据读数据

library(entropy)
gene_exp <- read.table("genefile.txt",sep=" ",header=T)
gene_exp <- as.matrix(gene_exp)
gene_id <- read.table("geneid.txt",header=T)

数据观察

表达数据是54675*22，id是54675与表达数据行数匹配

再看看课件里有啥别的信息

表达数据是11列加11列

大体思路

1. 表达稳定性（熵）如何得到

我是一个喜欢手写的原始人

2. 将各自的熵和gene_id对上号
略
3. 提取前10不稳定
按熵降序排序

代码和结果

1. 表达稳定性（熵）获取

giveMeYourEntropy <- function(x){
  n <- length(x)
  dis1 <- discretize(x[1:(n/2)],numBins = 3,r=range(x[1:(n/2)]))
  dis2 <- discretize(x[(n/2):n],numBins = 3,r=range(x[(n/2):n]))
  freqs2d <- rbind(dis1,dis2)
  en <- mi.plugin(freqs2d,unit="log2")
  return (en)
}

虽然我是已知11+11啦，但是写成长度一劈二去算两个离散化后频率其实是想体现我尊重代码！在函数里写一些已知的数字不会很没有美感嘛哈哈哈。

2. 将各自的熵和gene_id对上号

gene_processed <- gene_id
gene_processed$entropy <- apply(gene_exp,1,giveMeYourEntropy)

众所周知我是一条懒狗
没必要把熵的信息再加到表达那个大大大矩阵屁股上。既然id和表达数据的行一一对应，我为啥不直接把熵加载一列的id后面呢？

3. 提取前10不稳定

sorted <- gene_processed[order(gene_processed$entropy,decreasing=T),]
head(sorted,10)

前10不稳定