Создание графика QQ в R

Я пытался создать график qq в R. Мне было трудно использовать свои результаты, поэтому я попытался следовать примеру из «Базовый статистический анализ в генетических исследованиях случай-контроль, Clarke et al.)

Шаг 5, а, iii) Я заменил путь к моделям и обязательным полям и выглядит следующим образом:

data<-read.table("C:\Users\X\Desktop\BIOM3006\Alternate/data.assoc.logistic",header=TRUE);pdf("C:\Users\X\Desktop\BIOM3006\Alternate/pvalue.qq.plot.pdf");

obs<-−log10(sort(data[data$TEST=="ADD",]$P));exp<-−log10(c(1:length(obs))/(length(obs)+ 1));plot(exp, 

obs<-ylab="Observed(−logP)",xlab="Expected(−logP)",ylim=c(0,20),xlim=c(0,7)) lines(c(0,7),c(0,7),col=1,lwd=2);dev.off()

Это сообщение об ошибке, которое я получаю:

data<-read.table("C:\Users\Tom\Desktop\BIOM3006\Alternate/data.assoc.logistic",header=TRUE);pdf("C:\Users\Tom\Desktop\BIOM3006\Alternate/pvalue.qq.plot.pdf");
Error: '\U' used without hex digits in character string starting ""C:\U"

obs<-−log10(sort(data[data$TEST=="ADD",]$P));exp<-−log10(c(1:length(obs))/(length(obs)+ 1));plot(exp, 
Error in log10(sort(data[data$TEST == "ADD", ]$P)) : 
 non-numeric argument to mathematical function

In addition: Warning message:
In is.na(x) : is.na() applied to non-(list or vector) of type 'NULL'

obs<-ylab="Observed(−logP)",xlab="Expected(−logP)",ylim=c(0,20),xlim=c(0,7))
Error: unexpected ',' in "obs<-ylab="Observed(-logP)","

lines(c(0,7),c(0,7),col=1,lwd=2)
Error in plot.xy(xy.coords(x, y), type = type, ...) : 
plot.new has not been called yet

;dev.off()
Error: unexpected ';' in ";"

Я все еще разбираюсь с этим программным обеспечением, поэтому любая помощь будет оценена и извините, если я упустил из виду что-то явно очевидное. Том


r quantile genetics
person user3253698    schedule 05.04.2014    source источник
comment
Я бы сказал, перепроверьте путь к файлу, который вы используете для загрузки данных - похоже, у вас есть косые черты, идущие в обоих направлениях.   -  person Stedy    schedule 05.04.2014
comment
Стеди, наверное, прав. Используйте вывод .Platform$file.sep, чтобы увидеть, как разделить элементы вашего пути.   -  person alexis_laz    schedule 06.04.2014
comment
Это не имеет ничего общего с графиками Q-Q. Обратная косая черта в R используется для экранирования следующего символа. \U начинает определение символа Юникода, поэтому он ожидает шестнадцатеричное значение после \U. Попробуйте ввести, например, \U03A3 в консоли R. В путях к файлам используется косая черта (также и в Windows). См. этот пост.   -  person jlhoward    schedule 06.04.2014

Ответы (3)


arrow_upward
0
arrow_downward

Я бы сказал, что вы довольно близки, были некоторые проблемы с форматированием, которые требуют незначительных исправлений (например, отсутствие пробелов между <- и -log), но в остальном все выглядит нормально, хотя трудно сказать без воспроизводимого примера.

data <- read.table("C:\Users\X\Desktop\BIOM3006\Alternate/data.assoc.logistic",header=TRUE)
pdf("C:\Users\X\Desktop\BIOM3006\Alternate/pvalue.qq.plot.pdf")
obs <- −log10(sort(data[data$TEST=="ADD",]$P))
exp<- −log10(c(1:length(obs))/(length(obs)+ 1))
plot(exp, obs, ylab="Observed(−logP)",xlab="Expected(−logP)",ylim=c(0,20),xlim=c(0,7)) lines(c(0,7),c(0,7),col=1,lwd=2))
dev.off()
person Stedy    schedule 05.04.2014

arrow_upward
0
arrow_downward

Соглашусь с комментариями выше, проблема не в самом вашем куске кода, а в файле "\U".

Кстати, я публикую более простой скрипт для создания графиков QQ с выводом plink, как в вашем примере.

# If you call the script from bash it may include the name of your plink.qassoc file
args=(commandArgs(TRUE))
ARG = args[1]

# check libraries. It may download some from the Internet
if(!("methods" %in% rownames(installed.packages())))  
install.packages("methods")
library("methods")
if(!("ggplot2" %in% rownames(installed.packages())))  
install.packages("ggplot2")
library("ggplot2")

# Input data file (it is different for your file)
filename = paste ("../plink_files/", ARG, ".plink.qassoc.adjusted",sep="")
data = read.table (file = filename, header=T)

# Takes proability vector. 
PROB = data$UNADJ

# minus Log-P
logQQ = -log(data$QQ  ,10)
logP  = -log(PROB,10)

# QQ Plot
qqplot = ggplot (data = data, aes (x = logQQ, y = logP )) + geom_point(shape=1)
  qqplot = qqplot + xlab("Expected -logP values") + ylab("Observed -logP values")
  qqplot = qqplot + geom_abline (intercept = 0, slope = 1, color="red")
  qqplot = qqplot + ylim(0,max(logQQ,logP)) + xlim (0,max(logQQ,logP))

# Saving plot into a file
output = paste ("../results/", ARG, "_QQplot.jpg", sep="")
ggsave(output, plot = qqplot,width=6.47,height=4.67)
person elcortegano    schedule 19.04.2016

arrow_upward
0
arrow_downward

ДЛЯ МЕНЯ ЭТО РАБОТАЕТ: я просто делаю небольшие модификации. ПОСМОТРЕТЬ данные‹-read.table("data.assoc",header=TRUE); pdf("pvalue.qq.plot.pdf");

obs‹--log10(sort(data[data$TEST=="ADD",]$P)); exp‹--log10(c(1:длина(наблюдения))/(длина(наблюдения)+ 1));

график (exp, ylab = «Наблюдается (-logP)», xlab = «Ожидается (-logP)», ylim = c (0,20), xlim = c (0,7))

строки (с (0,7), с (0,7), столбец = 1, lwd = 2);

dev.off()

person forever    schedule 14.11.2016