First, I want to thank you for making a detailed and reproducible report.
 I have only done enough to exonerate the current version.

Second, there was a problem in the script you gave -- I see no way for
chrnames="chr1" to work with current annotation specifications as the
chrnames is used to query the geneannopk package.  for illuminaHumanv1.db
you have to use "1".

Third, the following additional code shows that the trans scores retrieved
involve SNP-gene distances > 100000 bp, USING THE BIOCONDUCTOR ANNOTATION
PACKAGES SPECIFIED IN THE CALL.

> tt = transTab(trS)  # get a data frame of the scores and some location
information
> dim(tt)
[1] 2853880       8

# that gives 20 per-SNP association scores, for the 20 largest trans
associations on chr1 to SNP on chr1, that is, the best among SNP-gene pairs
satisfying
# your distance restriction

> library(SNPlocs.Hsapiens.dbSNP.20111119)
> sl1 = getSNPlocs("ch1", as.GRanges=TRUE)
> rsn1 = paste("rs", values(sl1)$RefSNP_id, sep="")
>
> sst = start(sl1)
> names(sst) = rsn1
>
> MM <- match(as.character(tt[,1]), names(sst), nomatch=0)
> sst2 =sst[MM]
>
> tt2 = tt[which(MM>0),]
>
> summary(abs(abs(tt2$geneloc)-sst2))
     Min.   1st Qu.    Median      Mean   3rd Qu.      Max.      NA's
   100100  35100000  76200000  88130000 137300000 248300000     53284

this shows that the SNP-gene pairs on which scores were returned are
separated by at least 100100 b

it also shows that some gene/snp pairs couldn't be related by location, and
that needs to be rationalized.  it may be something
simple.

I think the problem that you report comes from using NCBI2R SNP annotation
to characterize SNP
locations ... I'm not familiar with that package but I assume that it
yields SNP locations not completely consistent with the
SNPlocs packaged used in the filltering for association testing.

It would be possible to use the NCBI2R locations for filtering if you made
GRanges with the names and locations corresponding
to the genotyped variants and put them in a package with a getSNPlocs
function that would emulate the SNPlocs
package behaviors.

Another approach would work from passing the "cis map" in as a list,
indicating which SNPs are excluded from
testing for each gene.  If these become compelling use cases I can
introduce this capability.

Let me know if there are further issues.

On Thu, Oct 4, 2012 at 12:09 PM, Vincent Carey
<stvjc@channing.harvard.edu>wrote:

> Thank you for this report.  I will get back to you ASAP.  This is a tough
> time with versions changing, but I wonder if you could try with the new
> release.   The robustification of location handling has been a key concern
> of development in the past couple of months.
>
> On Thu, Oct 4, 2012 at 11:56 AM, James Peters [guest] <
> guest@bioconductor.org> wrote:
>
>>
>> Hi there,
>>
>> I'm using the excellent GGtools to package to look for eQTLs. When I use
>> the transScore function to look for trans eQTL I find that some of the
>> eQTL-gene pairs I discover are closer than the distance specified by the
>> argument 'radius' (used to define what is cis; hence eSNPs closer to the
>> target gene than radius should not be returned by transScores).
>>
>> Any suggestions on what I'm doing wrong would be greatly appreciated.
>>
>> Best wishes
>>
>> James
>>
>> An example using GGdata looking only at genes and SNPs on chromosome 1.
>> Here I define trans as > 100,000 bps.
>>
>> library(GGtools)
>> library(GGdata)
>>
>> # get radius used to define cis
>> rad <- 1e+05
>>
>> trS <-transScores("GGdata", snpchr = as.character(1), rhs = ~1, K = 20,
>> targdirpref = "tsco", geneApply = lapply, chrnames = paste("chr",
>> as.character(1), sep = ""), radius = rad, renameChrs = NULL, probesToKeep =
>> NULL,
>> batchsize = 200, genegran = 50, shortfac = 10, wrapperEndo= function(x)
>> MAFfilter(x, low=0.1),
>> geneannopk = "illuminaHumanv1.db", snpannopk =
>> "SNPlocs.Hsapiens.dbSNP.20111119", gchrpref = "",
>> schrpref = "ch")
>>
>> # get gene indexes
>> TOPind <- geneIndcol(trS, 1)
>>
>> # construct results table using probe ids, chi squared #association
>> scores and snp loci
>> TOP <- data.frame(cbind(probeid= featureNames(getSS("GGdata",1))[TOPind],
>> score= topScores(trS), snpids = locusNames(trS)))
>>
>> TOP$score <- as.numeric(as.character(TOP$score))
>>
>> # set chi squared score threshold for significance for trans eQTL
>> (equivalent to p val < 1e-11)
>> threshold <- 46.328476
>>
>> tophits <- TOP[which(TOP$score > threshold),]
>>
>> # map probe to gene symbol, entrez id and chromosome
>> PROBE2SYMS <- select(illuminaHumanv1.db, keys =
>> as.character(tophits$probeid), cols = c("SYMBOL","ENTREZID","CHR"))
>>
>> tophits <- cbind(PROBE2SYMS, tophits)
>>
>> # get snp details
>> library(NCBI2R)
>> mysnps <- as.character(tophits$snpids)
>> annotSNPs <- AnnotateSNPList(mysnps, file="")
>>
>> selectSNPinfo <- annotSNPs[,c("marker","chr","chrpos")]
>>
>> # get annotation on the genes using NCBI
>> geneinfo <- GetGeneInfo(tophits$ENTREZID)
>> # ori gives info on whether the gene is on sense or #anti-sense strand
>> geneinfo <- geneinfo[,c("locusID","GeneLowPoint","GeneHighPoint","ori")]
>>
>> y <- cbind(tophits, selectSNPinfo)
>> #check
>> if(identical(y$marker, as.character(y$snpids))==FALSE){stop("The SNP ids
>> don't match up")}
>>
>> # merge gene info with snp info and trans eQTL results
>> final <- cbind(geneinfo, y)
>> # check
>> if(identical(final$locusID, final$ENTREZID)==FALSE){stop("The Entrez ids
>> don't match up")}
>>
>> # clean up the dataframe (remove duplicate cols, make col names clearer)
>> colnames(final)[1] <- 'EntrezID'
>> colnames(final)[8] <- 'geneCHR'
>> colnames(final)[13] <- 'snpCHR'
>> colnames(final)[14] <- 'snpPos'
>> final <- final[,-c(7,9,12)]
>>
>> # calculate distance from snp to gene start position
>> diff <- rep(NA, nrow(final))
>> for (i in 1:nrow(final)){
>>   if(final$ori[i] == '+'){
>>     diff[i] <- abs(final$GeneLowPoint[i] - final$snpPos[i])
>>   } else if
>>   (final$ori[i] == '-'){
>>     diff[i] <- abs(final$GeneHighPoint[i] - final$snpPos[i])
>>   }
>> }
>>
>> # eQTL which are in fact cis
>> final[which(diff < rad),]
>>
>>
>>
>>  -- output of sessionInfo():
>>
>> > sessionInfo()
>> R version 2.15.1 (2012-06-22)
>> Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit)
>>
>> locale:
>> [1] en_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8
>>
>> attached base packages:
>> [1] splines   stats4    stats     graphics  grDevices utils     datasets
>>  methods   base
>>
>> other attached packages:
>>  [1] SNPlocs.Hsapiens.dbSNP.20111119_0.99.6
>> SNPlocs.Hsapiens.dbSNP.20110815_0.99.6 GGdata_1.0.19
>>    illuminaHumanv1.db_1.12.2              org.Hs.eg.db_2.7.1
>>  [6] RSQLite_0.11.2                         DBI_0.2-5
>>          AnnotationDbi_1.18.4                   Biobase_2.16.0
>>             GGtools_4.4.0
>> [11] Rsamtools_1.8.6                        Biostrings_2.24.1
>>          GenomicRanges_1.8.13                   IRanges_1.14.4
>>             BiocGenerics_0.2.0
>> [16] GGBase_3.18.0                          snpStats_1.6.0
>>           Matrix_1.0-9                           lattice_0.20-10
>>              survival_2.36-14
>> [21] NCBI2R_1.4.3
>>
>> loaded via a namespace (and not attached):
>>  [1] annotate_1.34.1          biomaRt_2.12.0           bit_1.1-8
>>        bitops_1.0-4.1           BSgenome_1.24.0          ff_2.2-7
>>       genefilter_1.38.0
>>  [8] GenomicFeatures_1.8.3    grid_2.15.1              RCurl_1.95-0.1.2
>>       rtracklayer_1.16.3       tools_2.15.1
>> VariantAnnotation_1.2.11 XML_3.95-0
>> [15] xtable_1.7-0
>>
>> --
>> Sent via the guest posting facility at bioconductor.org.
>>
>
>

	[[alternative HTML version deleted]]

