劉小澤寫于18.12.10
生信必備三大件：生物斧蜕、統(tǒng)計批销、技術(shù)均芽，我想要借助R來學(xué)習(xí)統(tǒng)計學(xué)知識，因為平時使用R比較頻繁深纲，理解起來應(yīng)該也會更快一些湃鹊。先來一次R的知識迭代币呵，更新下知識庫
先看一些比較基礎(chǔ)的R操作
這其中有一些是我之前沒有學(xué)習(xí)到的余赢，第一次get感覺超級有用
我們一般都會查詢陌生的函數(shù)如何使用哈垢，會用到
help()妻柒。默認狀態(tài)下，help()函數(shù)只會在載入加載過的程序包中搜索耘分，如果想在所有包中進行搜索举塔，又不想先加載它，可以用help("xx", try.all.packages=TRUE)陶贼；如果已經(jīng)知道函數(shù)在哪個包中啤贩，可以用help("xx", package="yy")-
數(shù)據(jù)直接編輯[這個東西很有意思]：借助
edit()函數(shù)，我們可以直接在R中進行數(shù)據(jù)修改（這個功能我其實一直在尋找拜秧，但最近才發(fā)現(xiàn)，這樣就不用導(dǎo)出后修改然后再導(dǎo)入R了）枉氮。
例如志衍，我們要修改mtcars這個數(shù)據(jù)集，但是為了避免修改原有的數(shù)據(jù)集聊替，可以新保存到一個變量如new <- edit(mtcars)楼肪，然后修改完點quit就好；
如果想在現(xiàn)有數(shù)據(jù)集基礎(chǔ)上修改惹悄，不想新建春叫，那么可以用fix（mtcars）谢肾；
創(chuàng)建一個新的數(shù)據(jù)集舱权，再向其中填充數(shù)據(jù)new <- edit(data.frame())，默認創(chuàng)建兩列，然后拖動右下角就出現(xiàn)更多的列笼痹；Windows是可以直接使用edit的韭赘，Mac的話需要先安裝XQuartz這個工具
如果想研究數(shù)據(jù)集中的某個變量扁掸，比如mtcars中的mpg乓诽，可以使用
mtcars$mpg，但是如果還需要別的變量晨横，那還是多敲幾遍mtcars$洋腮；其實可以先將數(shù)據(jù)集激活，放到環(huán)境中手形，什么時候想提取啥供，就直接輸變量名就好了，比如先attach(mtcars)叁幢，就是將mtcars當(dāng)做當(dāng)前數(shù)據(jù)集滤灯，然后想調(diào)用mpg變量直接輸入mpg就好坪稽；不想用了就detach(mtcars)-
不同數(shù)據(jù)曼玩，不同目的，不同分析：
屬性數(shù)據(jù)
例如mtcars中的cyl數(shù)據(jù)（記錄氣缸數(shù)）是屬性變量窒百。先激活mtcars數(shù)據(jù)集黍判，然后使用
table(cyl)得到cyl中的三個值：4，6篙梢，8和相應(yīng)的頻數(shù)顷帖；還可以barplot畫直方圖barplot(table(cyl))，直接barplot(cyl)也能畫渤滞，但是分組信息就體現(xiàn)不出來贬墩，失去了意義數(shù)值型數(shù)據(jù)
這一類在統(tǒng)計分析中經(jīng)常用到，例如mtcars中的mpg變量就是數(shù)值型數(shù)據(jù)妄呕；另外不同函數(shù)反映的意義也是不同的陶舞。
作圖可以畫直方圖
hist(mpg)，箱線圖boxplot()绪励；統(tǒng)計分析可以分析均值
mean()肿孵；計算截掉5%的均值mean(mpg, trim=.05)；按分組變量cyl來計算mpg的分組均值tapply(mpg, cyl, mean)疏魏；計算cyl為6的mpg的均值mean(mpg[cyl == 6 ])停做；計算常用分位數(shù)（包括極小、極大大莫、中位數(shù)蛉腌、兩個四分位數(shù)）quantile(mpg)或者fivenum(mpg)或者summary(mpg)；計算四分位數(shù)的極差IQR(mpg)；標準差sd(mpg)烙丛；中位絕對離差(median absolute deviation)mad(mpg)探索二元關(guān)系
1 擬合線性回歸贯吓，例如
> z <- lm(cyl~mpg) > z # 結(jié)果 Call: lm(formula = cyl ~ mpg) Coefficients: (Intercept) mpg 11.2607 -0.2525 # 線性回歸截距11.26，斜率-0.252 相關(guān)系數(shù)：考察回歸擬合的好壞蜀变，例如相關(guān)系數(shù)R可以用
cor(cyl, mpg)表示悄谐，更常用的R2 就能計算出來為0.726，表示數(shù)據(jù)變化的72.6%可以用氣缸數(shù)（cyl）與每加侖的英里數(shù)（mpg）表示3 殘差分析：殘差是估計值與觀測值的偏差
z <- lm(cyl ~ mpg) #將回歸分析的結(jié)果作為對象保存在lm.res中 lm.resids <- resid(z) # 提取殘差向量 plot(lm.resids) # 畫個散點圖 hist(lm.resids) # 畫個直方圖库北，是否為鐘形爬舰？ qqnorm(lm.resids) # QQ圖(Quantile-Quantile Plots)是否落在直線上？ # 都是為了檢驗數(shù)據(jù)集分布假設(shè)是否有效
了解下R的對象
對象屬性
R運行都是靠對象寒瓦，所有的對象都有兩個內(nèi)在屬性：類型和長度
類型包括四種：數(shù)值型情屹、字符型、復(fù)數(shù)型杂腰、邏輯型（T/F/NA）垃你，用函數(shù)mode()查看。另外不管什么類型的數(shù)據(jù)喂很，缺失值都是用NA（Not Available）表示惜颇，不是數(shù)值用NaN(Not a Number)表示；
數(shù)值型：數(shù)值太大用指數(shù)形式表示少辣，如
N <- 2e20凌摄，Inf表示正無窮-
字符型：輸入時要加上雙引號，如果要在其中繼續(xù)引用雙引號的話漓帅，可以用
\進行轉(zhuǎn)義锨亏；或者直接用單引號> x <- "Double quotes \" delimitate R's strings." # 雙引號 > x <- 'Double quotes " delimitate R\'s strings.' # 單引號 >x [1] "Double quotes \" delimitate R's strings." # 在R中顯示的狀態(tài)[注意這里在R中顯示的還帶著轉(zhuǎn)義符，其實顯示出來是沒有的] > cat(x) Double quotes " delimitate R's strings. # 實際上的狀態(tài)
長度是對象中元素的數(shù)目忙干，用length()查看
對象類別
只有數(shù)據(jù)框和列表支持多種對象并存器予；因子只有兩種類型
向量是一個變量的取值；因子是一個分類變量捐迫；數(shù)組是一個k維的數(shù)據(jù)表乾翔；矩陣是數(shù)組的特例【數(shù)組或者矩陣的所有元素都是一種】；數(shù)據(jù)框是一個或幾個向量/因子構(gòu)成弓乙，并且它們必須等長末融，可以是不同的類型；列表可以包含任何類型的對象（也可以列表包含列表）
瀏覽對象
利用ls()可以查看當(dāng)前在內(nèi)存中的對象暇韧，但是這個函數(shù)只列出了對象名勾习，并且是所有的；
想要查看名稱中含有某個指定字符（如x）的對象懈玻，可以指定pattern：ls(pat = "x");
想要看以某個字母開頭的對象巧婶，可以利用ls(pat = "^x") ;
如果想看所有對象的詳細信息呢？ls.str()
刪除對象
-
rm(x,y)刪除對象x和y -
rm(list=ls())刪除所有對象 -
rm(list=ls(pat="^x"))刪除所有x開頭的對象
了解下R的向量
建立向量
數(shù)值型向量
向量沒規(guī)律：
c()-
向量的規(guī)律比較簡單：
seq()或“:”> 1:10 # 1 2 3 4 5 6 7 8 910 > 1:10-1 # 0 1 2 3 4 5 6 7 8 9 > 1:(10-1) # 注意有括號與上面沒有的區(qū)別 # 1 2 3 4 5 6 7 8 9 > z <- seq(1,5,by=0.5) # 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 > z <- seq(1,10,length=11) # 1.0 1.9 2.8 3.7 4.6 5.5 6.4 7.3 8.2 9.1 10.0 -
向量的規(guī)律比較復(fù)雜：
rep()> z <- rep(2:5,2) # 2 3 4 5 2 3 4 5 > z <- rep(2:5,each=2) # 2 2 3 3 4 4 5 5 -
通過鍵盤輸入向量
scan()> z <- scan() 1: 1 2 3 4 5 6 7: # Read 6 items > z # 1 2 3 4 5 6
字符型向量
注意引號在輸入時應(yīng)該寫作：\" ；paste（）可以連接多個參數(shù)成為字符串艺栈，其中如果有數(shù)值英岭，那么數(shù)值會被強制轉(zhuǎn)為字符串；默認空格分割各個字符湿右，使用sep=自定義
邏輯型向量
TRUE诅妹、FALSE可以簡寫成T、F毅人；如果轉(zhuǎn)換為數(shù)值吭狡，F(xiàn)ALSE為0，TRUE為1
因子型向量
利用factor（）建立
factor(x, levels = sort(unique(x), na.last = TRUE),
labels = levels, exclude = NA, ordered = is.ordered(x))
levels 指定因子水平（缺省值是向量x中不同的值）丈莺；
labels 指定水平名稱划煮；
ordered是邏輯型選項，表示因子水平是否有順序
- 字符型因子轉(zhuǎn)為數(shù)值型因子【反之亦然】
> a <- c("x","y","z","x")
> a
[1] "x" "y" "z" "x"
> levels(a) <- c(1,2,3)
> a
[1] 1 2 3 1
Levels: 1 2 3
-
產(chǎn)生規(guī)則的因子序列缔俄，利用
gl(k,n)弛秋，k是水平數(shù)，n是重復(fù)數(shù)俐载；另外還有兩個選項蟹略，length指定總共數(shù)據(jù)個數(shù)，label指定每個水平名稱> gl(2,3) [1] 1 1 1 2 2 2 Levels: 1 2 > gl(2,3,length=10) [1] 1 1 1 2 2 2 1 1 1 2 Levels: 1 2 > gl(2,3,labels = c("x","y")) [1] x x x y y y Levels: x y
向量提取
這里簡單說下根據(jù)邏輯值提取
# 例如x向量是這樣的
> x <- c(20, 15, 11, 6)
> x>10 # 這個僅僅是判斷瞎疼，返回的也是F/T
[1] TRUE TRUE TRUE FALSE
> x[x>10] # 提取大于10的元素
[1] 20 15 11
> x[ x < 15 & x > 10] # 限定多個提取條件
[1] 11
> y <- x[!is.na(x)] # 找到x中的非缺失值
> z <- x[(!is.na(x))&(x>0)] # 非負非缺失值
了解下R的矩陣 Matrix
相比于數(shù)組科乎，矩陣使用頻率更高，構(gòu)建矩陣使用matrix
> x <- matrix(1:3, 2, 3) # 默認按列填充
> x
[,1] [,2] [,3]
[1,] 1 3 2
[2,] 2 1 3
# 如果改成按行填充贼急，將 byrow=TRUE就好
矩陣元素的替換
> x[,3] <- NA
> x
[,1] [,2] [,3]
[1,] 1 3 NA
[2,] 2 1 NA
> x[is.na(x)] <- 1 # 缺失值用1代替
> x
[,1] [,2] [,3]
[1,] 1 3 1
[2,] 2 1 1
矩陣統(tǒng)計運算
矩陣的排列是有方向性的，規(guī)定矩陣按列排列捏萍，一般不說明的時候太抓，統(tǒng)計函數(shù)也是按列計算（但是可以用MARGIN來改變，等于1代表列令杈，等于2代表行）
cov()與cor()分別計算矩陣的協(xié)方差矩陣和相關(guān)系數(shù)陣走敌；
可以進行標準化scale(x, center=T, scale=T) ；
按列求均值apply(x, MARGIN=2, FUN=mean)
了解下R的數(shù)據(jù)框 Data frame
雖然說數(shù)據(jù)框與矩陣很相似逗噩，也是二維表格掉丽，也是要求各個變量的觀測值長度相等，但是异雁，在數(shù)據(jù)框中捶障，行和列的意義是不一樣的， 其中列表示變量纲刀，行為觀測
提取
- 提取一列：
數(shù)據(jù)框名$變量名 - 提取滿足條件子集：例如
subset(df, con1 == "treated" & con2 > 100)
添加
- 基本方法：
df$val1 <- condf是數(shù)據(jù)框项炼，val1是新變量，con是定義的val1元素 - 使用with：
df$val1 <- with(df, con) - transform一次增加多個變量：
df <- transform(df, val1 = con, val2 = con2)
數(shù)據(jù)保存
一般用write.table()或者save()
保存為簡單的文本文件：
write.table(df, file="/path/", row.names = F, quote = F)
row.names =F 表示行名不寫入文件；quote = F表示變量名不放在雙引號中保存為逗號分隔（csv）：用
write.csv()-
保存為R格式文件：
save(df, file="/path/")或者保存全部
save.image()锭部，它等價于save(list =ls(all=TRUE), file=".RData")
數(shù)據(jù)讀取
讀取數(shù)據(jù)集
R內(nèi)置的基本數(shù)據(jù)集有100多個（常為數(shù)據(jù)框和列表）暂论。它們隨R的啟動全部一次性自動載入，通過命令data() 可以查看全部的數(shù)據(jù)集（也包含了通過library()加載的包中數(shù)據(jù)集）拌禾；使用data(package = "pkname") 可以列出包pkname中的所有數(shù)據(jù)集取胎，但是可能還未被加載，確定要用的時候可以加載包
讀取Rdata
涉及多個數(shù)據(jù)集的分析時湃窍，最常使用load("/path/")
R的編程思維
控制結(jié)構(gòu)中扼菠，條件語句少不了，例如
if (條件) 表達式1 else 表達式2坝咐，但是這樣有些冗余循榆，可以利用ifelse(條件， 表達式1墨坚， 表達式2)循環(huán)結(jié)構(gòu)中秧饮，知道終止條件用for（）；不知道運行次數(shù)用while（）
-
一般將一組命令放在大括號中=》模塊化編程
# 一般格式：函數(shù)名 = function(變量列表) 函數(shù)體 # 例如一個畫圖函數(shù) > myfun <- function(x, y) { data <- read.table(x) plot(data$V1, data$V2, type="l") title(y) } 代碼要增加注釋泽篮，設(shè)置行前自動縮進
-
命令向量化（vectorization） 盗尸，就是將循環(huán)隱含在表達式中，例如條件語句可以用邏輯索引向量代替
for (i in 1:length(x)){ if (x[i] == b) x[i] <- 0 else y[i] <- 1 } # 這么長的函數(shù)帽撑，其實就是為了判斷x中元素是否為b泼各，等于b就賦值0；不等就賦值1 # 可以用下面代替 x[x == b] <- 0 x[x != b] <- 1向量化快的原因是：R中使用向量化會立即調(diào)用C語言亏拉，計算時比R快100倍
創(chuàng)建project來運行項目扣蜻，這樣一是可以方便讀取存儲；二是不同項目分類存放及塘，方便查找更新
想要把腳本從頭到尾運行一遍莽使，用
source()，記得添加絕對路徑
歡迎關(guān)注我們的公眾號～_～
我們是兩個農(nóng)轉(zhuǎn)生信的小碩笙僚，打造生信星球芳肌，想讓它成為一個不拽術(shù)語、通俗易懂的生信知識平臺肋层。需要幫助或提出意見請后臺留言或發(fā)送郵件到Bioplanet520@outlook.com
