RubyWoodward commented 1 year ago

Hello,

I am planning to on running bivariate rg analysis across all given 2495 genomic loci. This will entail testing ~735 traits against 2 phenotypes in the long run, although I am using a smaller number currently to test the method.

I have attempted to adapt the example analysis script and execute this in a bash script for just one trait against these two phenotypes, however, I get the error: "Error: Invalid target phenotype specified: 'OA'", which is subsequently repeated for my other phenotype.

I have ran several iterations of this, all resulting in a similar error. I understand that I may only enter one target phenotype at a time, which also results in an error within the function run.uni.bivar.

My adapted R script is as follows:

arg = commandArgs(T); ref.prefix = arg[1]; loc.file = arg[2]; info.file = arg[3]; target=arg[4]; out.fname = arg[5]

Load package

library(LAVA)

Read in data

loci = read.loci(loc.file); n.loc = nrow(loci) input = process.input(info.file, sample.overlap.file=NULL, ref.prefix)

Set univariate pvalue threshold

univ.p.thresh = 0.05

Analyse

print(paste("Starting LAVA analysis for",n.loc,"loci")) progress = ceiling(quantile(1:n.loc, seq(.05,1,.05))) # (if you want to print the progress)

u=b=list() for (i in 1:n.loc) { if (i %in% progress) print(paste("..",names(progress[which(progress==i)]))) # (printing progress) locus = process.locus(loci[i,], input)
if (!is.null(locus)) {

extract some general locus info for the output

            loc.info = data.frame(locus = locus$id, chr = locus$chr, start = locus$start, stop = locus$stop, n.snps = locus$n.snps, n.pcs = locus$K)  
            # run the univariate and bivariate tests
            loc.out = run.univ.bivar(locus, target, univ.thresh = univ.p.thresh)
            u[[i]] = cbind(loc.info, loc.out$univ)
            if(!is.null(loc.out$bivar)) b[[i]] = cbind(loc.info, loc.out$bivar)
    }

}

save the output

write.table(do.call(rbind,u), paste0(out.fname,".univ.lava"), row.names=F,quote=F,col.names=T) write.table(do.call(rbind,b), paste0(out.fname,".bivar.lava"), row.names=F,quote=F,col.names=T)

print(paste0("Done! Analysis output written to ",out.fname,".*.lava"))

I assume it is the run.uni.bivar function in which is causing these issues, I have attempted to edit this function within R itself to fit my situation, but no luck so far. I will be extremely grateful for any advice as to how I can advance further.

RubyWoodward commented 1 year ago

Hello, I have managed to work this out, maybe I missed something obvious!

zakathornton commented 9 months ago

Hi @RubyWoodward, I am having a similar issue to this - I get the error you got. How did you resolve it?

josefin-werme / LAVA