首页
学习
活动
专区
圈层
工具
发布
社区首页 >问答首页 >as.vector(x)中的错误:没有将此S4类胁迫到向量的方法

as.vector(x)中的错误:没有将此S4类胁迫到向量的方法
EN

Stack Overflow用户
提问于 2021-05-03 16:18:58
回答 1查看 3.4K关注 0票数 1

我试图在bash命令行上运行一个R脚本(我使用的是CentOS 8),命令是:cat 1_myScript.R | R --slave --args $SAMPLE"_x" $SAMPLE"_y"

其中$SAMPLE是我在R脚本中指定的参数,如下所示

> args<-commandArgs()

> aaa<-args[4]

这个语法总是适用于我的所有脚本,但是现在它给了我以下错误:

as.vector(x)中的错误:没有强制此S4类执行向量调用的方法: setdiff -> setdiff.default -> -> as.vector执行暂停

奇怪的是,如果我试图在R控制台中运行这个脚本

>source("1_myScript.R")

它没有任何错误。我已经查过了,它似乎是链接到我在脚本中使用的库"GenomicRanges“的一个特性。下面是我的脚本的正文(请注意,我不知道它失败的确切行):

代码语言:javascript
复制
#!/usr/bin/env Rscript

library(rtracklayer)
library(data.table)
library(tidyverse)
#args <- commandArgs()
#aaa<-args[4]
aaa<-gsub("_T.finalSorted.bam_CNVs","",args[1])
targDir<-"/srv/ngsdata/dalteriog/SV_analysis/example/WGS_NB_novogene/CNVDir/"
dataTable <-fread((paste0(targDir,args[2]), header=TRUE)
ratio<-data.frame(dataTable)

dataTable <-fread(paste0(targDir,args[1]), header=FALSE)
cnvs<- data.frame(dataTable)

ratio$Ratio[which(ratio$Ratio==-1)]=NA

cnvs.bed=GRanges(cnvs[,1],IRanges(cnvs[,2],cnvs[,3])) # primo ogetto del GRange obj --> chr; secondo --> range
ratio.bed=GRanges(ratio$Chromosome,IRanges(ratio$Start,ratio$Start),score=ratio$Ratio) # score è un metadata

overlaps <- subsetByOverlaps(ratio.bed,cnvs.bed) # regioni overlappanti i due df
normals <- setdiff(ratio.bed,cnvs.bed) # regioni diverse
normals <- subsetByOverlaps(ratio.bed,normals) # la stessa cosa, ma con lo score associato

#mu <- mean(score(normals),na.rm=TRUE)
#sigma<- sd(score(normals),na.rm=TRUE)

#hist(score(normals),n=500,xlim=c(0,2))
#hist(log(score(normals)),n=500,xlim=c(-1,1))

#shapiro.test(score(normals)[which(!is.na(score(normals)))][5001:10000])
#qqnorm (score(normals)[which(!is.na(score(normals)))],ylim=(c(0,10)))
#qqline(score(normals)[which(!is.na(score(normals)))], col = 2)

#shapiro.test(log(score(normals))[which(!is.na(score(normals)))][5001:10000])
#qqnorm (log(score(normals))[which(!is.na(score(normals)))],ylim=(c(-6,10)))
#qqline(log(score(normals))[which(!is.na(score(normals)))], col = 2)

numberOfCol=length(cnvs)

for (i in c(1:length(cnvs[,1]))) {
  values <- score(subsetByOverlaps(ratio.bed,cnvs.bed[i])) #score bayesiano della iesima CNV che overlappa con il file ratio 
  #wilcox.test(values,mu=mu)
  W <- function(values,normals){resultw <- try(wilcox.test(values,score(normals)), silent = TRUE)
    if(class(resultw)=="try-error") return(list("statistic"=NA,"parameter"=NA,"p.value"=NA,"null.value"=NA,"alternative"=NA,"method"=NA,"data.name"=NA)) else resultw}
  KS <- function(values,normals){resultks <- try(ks.test(values,score(normals)), silent = TRUE)
    if(class(resultks)=="try-error") return(list("statistic"=NA,"p.value"=NA,"alternative"=NA,"method"=NA,"data.name"=NA)) else resultks}
  #resultks <- try(KS <- ks.test(values,score(normals)), silent = TRUE)
  # if(class(resultks)=="try-error") NA) else resultks
  cnvs[i,numberOfCol+1]=W(values,normals)$p.value
  cnvs[i,numberOfCol+2]=KS(values,normals)$p.value
  }

if (numberOfCol==5) {
  names(cnvs)=c("chr","start","end","copy number","status","WilcoxonRankSumTestPvalue","KolmogorovSmirnovPvalue")  
}
if (numberOfCol==7) {
  names(cnvs)=c("chr","start","end","copy number","status","genotype","uncertainty","WilcoxonRankSumTestPvalue","KolmogorovSmirnovPvalue")  
}
if (numberOfCol==9) {
  names(cnvs)=c("chr","start","end","copy number","status","genotype","uncertainty","somatic/germline","precentageOfGermline","WilcoxonRankSumTestPvalue","KolmogorovSmirnovPvalue")  
}

cnvs$Wfdr <- p.adjust(cnvs$WilcoxonRankSumTestPvalue, method="BH",n=nrow(cnvs))
cnvs$KSfdr <- p.adjust(cnvs$KolmogorovSmirnovPvalue, method="BH",n=nrow(cnvs))
cnvs<-subset(cnvs, Wfdr <= 0.05 & KSfdr <= 0.05)



samp<-gsub("_T.finalSorted.bam_CNVs","",aaa)
cnvs<-add_column(cnvs, Sample=samp, .before=1)
write.table(cnvs, file=paste(targDir,aaa,"CNV.p.filtered.txt",sep=""),sep="\t",quote=F,row.names=F)
EN

回答 1

Stack Overflow用户

回答已采纳

发布于 2021-05-04 07:48:37

看来问题就在第22行:

代码语言:javascript
复制
normals <- setdiff(ratio.bed,cnvs.bed)

在本例中,我们希望找到两个granges对象之间的差异,但是执行此操作的函数是setdiff of GenomicRanges,它被包基所掩盖。答案的解决方案是这样编辑这样的行:

代码语言:javascript
复制
normals <- GenomicRanges::setdiff(ratio.bed,cnvs.bed)

再次感谢您的努力!

票数 1
EN
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接:

https://stackoverflow.com/questions/67372347

复制
相关文章

相似问题

领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档