bb45615b2f096a51c523f1b74f323893f054d0d8 angie Fri Jul 25 11:49:13 2014 -0700 For SNPs tracks option to append alleles to the rs# ID in the item name, by default show alleles on the + strand of reference.If a user still wants to see alleles on dbSNP's arbitrary and unstable strand, there's a new option to do so. refs #13694 diff --git src/hg/hgTracks/variation.c src/hg/hgTracks/variation.c index d6229f1..4b85511 100644 --- src/hg/hgTracks/variation.c +++ src/hg/hgTracks/variation.c @@ -444,120 +444,163 @@ struct bed4 *b = *((struct bed4 **)vb); int dif; if(a==0||b==0) return 0; dif = differentWord(a->chrom, b->chrom); if (dif == 0) dif = a->chromStart - b->chromStart; if (dif == 0) dif = a->chromEnd - b->chromEnd; if (dif == 0) dif = differentWord(a->name, b->name); return dif; } -static char *snp125ExtendName(char *rsId, char *chimpAllele, char *observedAlleles) +static char *snp125ExtendName(char *rsId, char *chimpAllele, char *observedAlleles, boolean doRc) /* Allocate and return a string of the format "rsId chimpAllele>observedAlleles" if - * chimpAllele is non-empty, otherwise "rsId observedAlleles". */ + * chimpAllele is non-empty, otherwise "rsId observedAlleles". If doRc, reverse-complement + * all alleles. */ { struct dyString *dy = dyStringNew(64); + +dyStringPrintf(dy, "%s ", rsId); if (isNotEmpty(chimpAllele)) - dyStringPrintf(dy, "%s %s>%s", rsId, chimpAllele, observedAlleles); + { + int chimpLen = strlen(chimpAllele); + if (doRc && isAllNt(chimpAllele, chimpLen)) + { + char chimpCopy[chimpLen+1]; + safecpy(chimpCopy, sizeof(chimpCopy), chimpAllele); + reverseComplement(chimpCopy, chimpLen); + dyStringPrintf(dy, "%s>", chimpCopy); + } else - dyStringPrintf(dy, "%s %s", rsId, observedAlleles); + dyStringPrintf(dy, "%s>", chimpAllele); + } +if (doRc) + { + int obsLen = strlen(observedAlleles); + char obsCopy[obsLen+1]; + safecpy(obsCopy, sizeof(obsCopy), observedAlleles); + char *obsWords[obsLen]; + int obsCount = chopByChar(obsCopy, '/', obsWords, obsLen); + int i; + // dbSNP orders alleles alphabetically so reverse the order too: + for (i = obsCount-1; i >= 0; i--) + { + char *al = obsWords[i]; + int alLen = strlen(al); + if (isAllNt(al, alLen)) + reverseComplement(al, alLen); + if (i < obsCount-1) + dyStringAppendC(dy, '/'); + dyStringPrintf(dy, "%s", al); + } + } +else + dyStringPrintf(dy, "%s", observedAlleles); return dyStringCannibalize(&dy); } -static void setSnp125ExtendedNameObserved(struct snp125 *snpList) -/* Append observed alleles to each snp's name. */ +static void setSnp125ExtendedNameObserved(struct snp125 *snpList, boolean useRefStrand) +/* Append observed alleles to each snp's name, possibly reverse-complementing alleles + * that were reported on the - strand. */ { struct snp125 *snpItem = snpList; for (; snpItem != NULL; snpItem = snpItem->next) - snpItem->name = snp125ExtendName(snpItem->name, NULL, snpItem->observed); + { + boolean doRc = useRefStrand && snpItem->strand[0] == '-'; + snpItem->name = snp125ExtendName(snpItem->name, NULL, snpItem->observed, doRc); + } } void setSnp125ExtendedNameExtra(struct track *tg) /* add extra text to be drawn in snp name field. This works by walking through two sorted lists and updating the name value for the SNP list with data from a table of orthologous state information */ { char cartVar[512]; safef(cartVar, sizeof(cartVar), "%s.extendedNames", tg->tdb->track); boolean enabled = cartUsualBoolean(cart, cartVar, // Check old cart var name for backwards compatibility w/ old sessions: cartUsualBoolean(cart, "snp125ExtendedNames", FALSE)); if (!enabled) return; + +safef(cartVar, sizeof(cartVar), "%s.allelesDbSnpStrand", tg->tdb->track); +boolean useRefStrand = ! cartUsualBoolean(cart, cartVar, FALSE); + struct sqlConnection *conn = hAllocConn(database); int rowOffset = 0; char **row = NULL; struct orthoBed *orthoItemList = NULL; /* list of orthologous state info */ struct orthoBed *orthoItem = orthoItemList; char *orthoTable = snp125OrthoTable(tg->tdb, NULL); struct sqlResult *sr = NULL; int cmp = 0; /* if orthologous info is not available, show only observed alleles */ if(isEmpty(orthoTable) || !sqlTableExists(conn, orthoTable)) { - setSnp125ExtendedNameObserved((struct snp125 *)tg->items); + setSnp125ExtendedNameObserved((struct snp125 *)tg->items, useRefStrand); hFreeConn(&conn); return; } /* get list of orthologous alleles */ sr = hRangeQuery(conn, orthoTable, chromName, winStart, winEnd, NULL, &rowOffset); while ((row = sqlNextRow(sr)) != NULL) { orthoItem = orthoBedLoad(row + rowOffset); if (orthoItem) slAddHead(&orthoItemList, orthoItem); } /* List of SNPs is already sorted, so sort list of Ortho info */ slSort(&orthoItemList, snpOrthoCmp); /* Walk through two sorted lists together */ struct snp125 *snpItem = tg->items; orthoItem = orthoItemList; while (snpItem!=NULL && orthoItem!=NULL) { /* check to see that we're not at the end of either list and that * the two list elements represent the same human position */ cmp = snpOrthoCmp(&snpItem, &orthoItem); + boolean doRc = useRefStrand && snpItem->strand[0] == '-'; if (cmp < 0) { // Update snp->name with observed alleles even if we don't have ortho data - snpItem->name = snp125ExtendName(snpItem->name, NULL, snpItem->observed); + snpItem->name = snp125ExtendName(snpItem->name, NULL, snpItem->observed, doRc); snpItem = snpItem->next; continue; } if (cmp > 0) { orthoItem = orthoItem->next; continue; } /* update the snp->name with the ortho data */ - snpItem->name = snp125ExtendName(snpItem->name, orthoItem->chimp, snpItem->observed); + snpItem->name = snp125ExtendName(snpItem->name, orthoItem->chimp, snpItem->observed, doRc); /* increment the list pointers */ snpItem = snpItem->next; orthoItem = orthoItem->next; } // If orthoItemList ends before snpItemList, add observed alleles to remaining snps: if (snpItem != NULL) - setSnp125ExtendedNameObserved(snpItem); + setSnp125ExtendedNameObserved(snpItem, useRefStrand); sqlFreeResult(&sr); hFreeConn(&conn); } static char *snp125MapItemName(struct track *tg, void *item) /* Now that snp125->name is overwritten when adding chimp allele suffix, we need * to strip it back off for links to hgc. */ { static char mapName[256]; struct snp125 *snp = item; safecpy(mapName, sizeof(mapName), snp->name); char *ptr = strchr(mapName, ' '); if (ptr != NULL) *ptr = '\0'; return mapName;