Презентация http://makarich.fbb.msu.ru/artemov/R/R_3.pdf

Парный t-test

(то, что осталось с предыдущего занятия, предыдущая презентация http://makarich.fbb.msu.ru/artemov/R/R_2.pdf )

Загрузим таблицу: пациент - изменение длительности сна при приёме лекарства1 - изменение длительности сна при приёме лекарства2

Является ли одно из лекарств лучше другого?

sleep.paired = read.table("http://makarich.fbb.msu.ru/artemov/R/sleep.paired.tab")
boxplot(sleep.paired$Drug1, sleep.paired$Drug2, names=c("Drug1", "Drug2"), col=c("red", "blue"))

Cпособ 1: парный t-test

sleep.test = t.test(sleep.paired$Drug1, sleep.paired$Drug2, paired=T)

Способ 2: рассмотрим разность sleep.paired$Drug2 - sleep.paired$Drug1. t.test с одним вектором на вход проверяет, отличается ли среднее в этом векторе значимо от 0.

diff = sleep.paired$Drug2 - sleep.paired$Drug1
sleep.test = t.test(diff)

Тесты ассоциации и циклы

Исходные данные

Таблица variant_counts.csv http://makarich.fbb.msu.ru/artemov/R/variant_counts.csv

Аргумент row.names=1 сообщает, что первый столбец таблицы хранит названия (номера) рядов.

varcounts=read.csv("http://makarich.fbb.msu.ru/artemov/R/variant_counts.csv", header=T, row.names=1)

Тесты ассоциации

Создадим вручную таблицу (матрицу) сопряженности:

gwas=matrix(c(3324,1896,2676,2104), nrow=2, ncol=2)
colnames(gwas)=c("reference", "mutant")
rownames(gwas)=c("Healthy", "Diseased")
gwas

Применим тест ассоциации:

fisher.test(gwas)
chisq.test(gwas)

Как вытащить p-value (число):

y=fisher.test(gwas)
y$p.value

Циклы

Посчитаем для каждой мутации частоту в европейской популяции. Для программ длиннее 1 строчки удобнее сохранять их исходный ход в файле: File -> New file -> R script

С символа # начинается комментарий - то, что игнорируется интерпретатором

#создадим пустой вектор, куда будем класть результаты
freq=c()
#1:nrow(varcounts) - вектор от 1 до кол-ва строк в varcounts. 
#Переменная i пробегает все значения этого вектора по очереди, для каждого выполняется то, что написано в фигурных скобках
for(i in 1:nrow(varcounts)){
    frequency_i=varcounts[i, "EA_alt"] / (varcounts[i, "EA_ref"] + varcounts[i, "EA_alt"])
    #добавляем вычисленное число в вектор результатов:
    freq=c(freq, frequency_i)
}

Функции и apply

Возьмем "числовую" часть dataframe-а и превратим его в матрицу:

mvarcounts=as.matrix(varcounts[,3:6])
typeof(mvarcounts)

Определим функцию, которая, получая на вход строку таблицы (вектор из 4 чисел), выполняет тест ассоциации и возвращает p-value:

association.test <- function(v){
        m=matrix(v, nrow=2, ncol=2) #создаем матрицу
        res=fisher.test(m)
        return(res$p.value)
}

Что такое "возвращает":

association.test(c(3,1,1,3))

Применим данную функцию к каждой строке матрицы - получим вектор p-value:

Второй аргумент - измерение таблицы (1 - строки, 2 - столбцы)

pv=apply(mvarcounts, 1, association.test)
varcounts$pvalue=pv

Kodomo

Пользователь