Java Code Examples for htsjdk.samtools.SAMReadGroupRecord#setSample()

/** Creates a simple header with the values provided on the command line. */
public SAMFileHeader createSamFileHeader() {
    final SAMReadGroupRecord rgroup = new SAMReadGroupRecord(this.READ_GROUP_NAME);
    rgroup.setSample(this.SAMPLE_NAME);
    if (this.LIBRARY_NAME != null) rgroup.setLibrary(this.LIBRARY_NAME);
    if (this.PLATFORM != null) rgroup.setPlatform(this.PLATFORM);
    if (this.PLATFORM_UNIT != null) rgroup.setPlatformUnit(this.PLATFORM_UNIT);
    if (this.SEQUENCING_CENTER != null) rgroup.setSequencingCenter(SEQUENCING_CENTER);
    if (this.PREDICTED_INSERT_SIZE != null) rgroup.setPredictedMedianInsertSize(PREDICTED_INSERT_SIZE);
    if (this.DESCRIPTION != null) rgroup.setDescription(this.DESCRIPTION);
    if (this.RUN_DATE != null) rgroup.setRunDate(this.RUN_DATE);
    if (this.PLATFORM_MODEL != null) rgroup.setPlatformModel(this.PLATFORM_MODEL);
    if (this.PROGRAM_GROUP != null) rgroup.setProgramGroup(this.PROGRAM_GROUP);

    final SAMFileHeader header = new SAMFileHeader();
    header.addReadGroup(rgroup);

    for (final String comment : COMMENT) {
        header.addComment(comment);
    }

    header.setSortOrder(this.SORT_ORDER);
    return header ;
}
 

private FingerprintIdDetails createUnknownFP(final Path samFile, final SAMRecord rec) {
    final PicardException e = new PicardException("Found read with no readgroup: " + rec.getReadName() + " in file: " + samFile);
    if (validationStringency != ValidationStringency.STRICT) {
        final SAMReadGroupRecord readGroupRecord = new SAMReadGroupRecord("<UNKNOWN>:::" + samFile.toUri().toString());
        readGroupRecord.setLibrary("<UNKNOWN>");
        readGroupRecord.setSample(defaultSampleID);
        readGroupRecord.setPlatformUnit("<UNKNOWN>.0.ZZZ");

        if (validationStringency != ValidationStringency.SILENT && missingRGFiles.add(samFile)) {
            log.warn(e.getMessage());
            log.warn("further messages from this file will be suppressed");
        }

        return new FingerprintIdDetails(readGroupRecord, samFile.toUri().toString());
    } else {
        log.error(e.getMessage());
        throw e;
    }
}
 

/**
 * Create a new HaplotypeBAMDestination
 *
 * @param sourceHeader SAMFileHeader used to seed the output SAMFileHeader for this destination.
 * @param haplotypeReadGroupID read group ID used when writing haplotypes as reads
 */
protected HaplotypeBAMDestination(SAMFileHeader sourceHeader, final String haplotypeReadGroupID) {
    Utils.nonNull(sourceHeader, "sourceHeader cannot be null");
    Utils.nonNull(haplotypeReadGroupID, "haplotypeReadGroupID cannot be null");
    this.haplotypeReadGroupID = haplotypeReadGroupID;

    bamOutputHeader = new SAMFileHeader();
    bamOutputHeader.setSequenceDictionary(sourceHeader.getSequenceDictionary());
    bamOutputHeader.setSortOrder(SAMFileHeader.SortOrder.coordinate);

    final List<SAMReadGroupRecord> readGroups = new ArrayList<>();
    readGroups.addAll(sourceHeader.getReadGroups()); // include the original read groups

    // plus an artificial read group for the haplotypes
    final SAMReadGroupRecord rgRec = new SAMReadGroupRecord(getHaplotypeReadGroupID());
    rgRec.setSample(haplotypeSampleTag);
    rgRec.setSequencingCenter("BI");
    readGroups.add(rgRec);
    bamOutputHeader.setReadGroups(readGroups);
    final List<SAMProgramRecord> programRecords = new ArrayList<>(sourceHeader.getProgramRecords());
    programRecords.add(new SAMProgramRecord("HaplotypeBAMWriter"));
    bamOutputHeader.setProgramRecords(programRecords);
}
 

protected SAMReadGroupRecord createReadGroupRecord(
        String RGID, String RGLB, String RGPL, 
        String RGPU, String RGSM, String RGCN, 
        String RGDS, Iso8601Date RGDT, Integer RGPI) {
    SAMReadGroupRecord rg = new SAMReadGroupRecord(RGID);
    rg.setLibrary(RGLB);
    rg.setPlatform(RGPL);
    rg.setSample(RGSM);
    rg.setPlatformUnit(RGPU);
    if(RGCN != null)
        rg.setSequencingCenter(RGCN);
    if(RGDS != null)
        rg.setDescription(RGDS);
    if(RGDT != null)
        rg.setRunDate(RGDT);
    if(RGPI != null)
        rg.setPredictedMedianInsertSize(RGPI);
    return rg;
}
 

public List<SAMReadGroupRecord> getReadGroups(final List<String> cellBarcodes) {
	List<SAMReadGroupRecord> g = new ArrayList<>(cellBarcodes.size());
	for (String id: cellBarcodes) {
		SAMReadGroupRecord rg = new SAMReadGroupRecord(id);
		rg.setLibrary(id);
	    rg.setPlatform(id);
	    rg.setSample(id);
	    rg.setPlatformUnit(id);
		g.add(rg);
	}
	return (g);
}
 

@BeforeClass
public void setUp() {
    header = ArtificialReadUtils.createArtificialSamHeader(1, 1, 1000);
    rg = new SAMReadGroupRecord(RGID);
    rg.setSample(sample);
    header.addReadGroup(rg);
    parser = new GenomeLocParser(header.getSequenceDictionary());
}
 

public void setupTest1(final int ID, final String readGroupId, final SAMReadGroupRecord readGroupRecord, final String sample,
                  final String library, final SAMFileHeader header, final SAMRecordSetBuilder setBuilder)
        throws IOException {

    final String separator = ":";
    final int contig1 = 0;
    final int contig2 = 1;
    readGroupRecord.setSample(sample);
    readGroupRecord.setPlatform(platform);
    readGroupRecord.setLibrary(library);
    readGroupRecord.setPlatformUnit(readGroupId);
    header.addReadGroup(readGroupRecord);
    setBuilder.setReadGroup(readGroupRecord);
    setBuilder.setUseNmFlag(true);

    setBuilder.setHeader(header);

    final int max = 800;
    final int min = 1;
    final Random rg = new Random(5);

    //add records that align to chrM and O but not N
    for (int i = 0; i < NUM_READS; i++) {
        final int start = rg.nextInt(max) + min;
        final String newReadName = READ_NAME + separator + ID + separator + i;

        if (i != NUM_READS - 1) {
            setBuilder.addPair(newReadName, contig1, start + ID, start + ID + LENGTH);
        } else {
            setBuilder.addPair(newReadName, contig2, start + ID, start + ID + LENGTH);
        }
    }
}
 

public void setupTest2(final int ID, final String readGroupId, final SAMReadGroupRecord readGroupRecord, final String sample,
                       final String library, final SAMFileHeader header, final SAMRecordSetBuilder setBuilder)
        throws IOException {

    final String separator = ":";
    final int contig1 = 0;
    final int contig2 = 1;
    final int contig3 = 2;
    readGroupRecord.setSample(sample);
    readGroupRecord.setPlatform(platform);
    readGroupRecord.setLibrary(library);
    readGroupRecord.setPlatformUnit(readGroupId);
    setBuilder.setReadGroup(readGroupRecord);
    setBuilder.setUseNmFlag(true);

    setBuilder.setHeader(header);

    final int max = 800;
    final int min = 1;
    final Random rg = new Random(5);

    //add records that align to all 3 chr in reference file
    for (int i = 0; i < NUM_READS; i++) {
        final int start = rg.nextInt(max) + min;
        final String newReadName = READ_NAME + separator + ID + separator + i;

        if (i<=NUM_READS/3) {
            setBuilder.addPair(newReadName, contig1, start + ID, start + ID + LENGTH);
        } else if (i< (NUM_READS - (NUM_READS/3))) {
            setBuilder.addPair(newReadName, contig2, start + ID, start + ID + LENGTH);
        } else {
            setBuilder.addPair(newReadName, contig3, start + ID, start + ID + LENGTH);
        }
    }
}
 

void setup(final int numReads,
           final String readName,
           final int ID,
           final String readGroupId,
           final SAMReadGroupRecord readGroupRecord,
           final String sample,
           final String library,
           final SAMFileHeader header,
           final SAMRecordSetBuilder setBuilder) throws IOException {
    final String separator = ":";
    readGroupRecord.setSample(sample);
    readGroupRecord.setPlatform(platform);
    readGroupRecord.setLibrary(library);
    readGroupRecord.setPlatformUnit(readGroupId);
    header.addReadGroup(readGroupRecord);
    setBuilder.setReadGroup(readGroupRecord);
    setBuilder.setUseNmFlag(true);
    setBuilder.setHeader(header);

    final int max = 15000;
    final int min = 1;
    final Random rg = new Random(5);

    for (int i = 0; i < numReads; i++) {
        final int start = rg.nextInt(max) + min;
        final String newReadName = readName + separator + ID + separator + i;
        setBuilder.addPair(newReadName, 0, start + ID, start + ID + 99);
    }
}
 

@Test(enabled = true, dataProvider = "LIBS_NotHoldingTooManyReads")
public void testLIBS_NotHoldingTooManyReads(final int nReadsPerLocus, final int downsampleTo, final int payloadInBytes) {
    logger.warn(String.format("testLIBS_NotHoldingTooManyReads %d %d %d", nReadsPerLocus, downsampleTo, payloadInBytes));
    final int readLength = 10;

    final SAMFileHeader header = ArtificialReadUtils.createArtificialSamHeader(1, 1, 100000);
    final int nSamples = 1;
    final List<String> samples = new ArrayList<>(nSamples);
    for ( int i = 0; i < nSamples; i++ ) {
        final SAMReadGroupRecord rg = new SAMReadGroupRecord("rg" + i);
        final String sample = "sample" + i;
        samples.add(sample);
        rg.setSample(sample);
        rg.setPlatform(NGSPlatform.ILLUMINA.getDefaultPlatform());
        header.addReadGroup(rg);
    }

    final boolean downsample = downsampleTo != -1;
    final DownsamplingMethod downsampler = downsample
            ? new DownsamplingMethod(DownsampleType.BY_SAMPLE, downsampleTo, null)
            : new DownsamplingMethod(DownsampleType.NONE, null, null);

    final WeakReadTrackingIterator iterator = new WeakReadTrackingIterator(nReadsPerLocus, readLength, payloadInBytes, header);

    final LocusIteratorByState li;
    li = new LocusIteratorByState(
            iterator,
            downsampler,
            false,
            samples,
            header,
            true
    );

    while ( li.hasNext() ) {
        final AlignmentContext next = li.next();
        Assert.assertTrue(next.getBasePileup().size() <= downsampleTo, "Too many elements in pileup " + next);
        // TODO -- assert that there are <= X reads in memory after GC for some X
    }
}
 

@BeforeClass
public void setUp() {
    header = ArtificialReadUtils.createArtificialSamHeader(1, 1, 1000);
    rg = new SAMReadGroupRecord(RGID);
    rg.setSample(sample);
    header.addReadGroup(rg);
    parser = new GenomeLocParser(header.getSequenceDictionary());
}
 

public static SAMFileHeader createSamHeader(final String sampleName){
    final SAMFileHeader samHeader = ArtificialReadUtils.createArtificialSamHeader(
            DEFAULT_NUM_CHROMOSOMES, DEFAULT_STARTING_CHROM, DEFAULT_CHROMOSOME_SIZE);
    final SAMReadGroupRecord readGroupRecord = new SAMReadGroupRecord(DEFAULT_READ_GROUP_NAME);
    readGroupRecord.setSample(sampleName);
    samHeader.addReadGroup(readGroupRecord);
    return samHeader;
}
 

private static SAMFileHeader initSAMFileHeader() {
    final SAMFileHeader samHeader = createArtificialSamHeader();
    SAMReadGroupRecord readGroup = new SAMReadGroupRecord(readGroupName);
    readGroup.setSample(DEFAULT_SAMPLE_NAME);
    samHeader.addReadGroup(readGroup);
    return samHeader;
}
 

protected int doWork() {
    IOUtil.assertInputIsValid(INPUT);
    IOUtil.assertFileIsWritable(OUTPUT);

    final SamReader in = SamReaderFactory.makeDefault()
        .referenceSequence(REFERENCE_SEQUENCE)
        .open(SamInputResource.of(INPUT));

    // create the read-group we'll be using
    final SAMReadGroupRecord rg = new SAMReadGroupRecord(RGID);
    rg.setLibrary(RGLB);
    rg.setPlatform(RGPL);
    rg.setSample(RGSM);
    rg.setPlatformUnit(RGPU);
    if (RGCN != null) rg.setSequencingCenter(RGCN);
    if (RGDS != null) rg.setDescription(RGDS);
    if (RGDT != null) rg.setRunDate(RGDT);
    if (RGPI != null) rg.setPredictedMedianInsertSize(RGPI);
    if (RGPG != null) rg.setProgramGroup(RGPG);
    if (RGPM != null) rg.setPlatformModel(RGPM);
    if (RGKS != null) rg.setKeySequence(RGKS);
    if (RGFO != null) rg.setFlowOrder(RGFO);

    log.info(String.format("Created read-group ID=%s PL=%s LB=%s SM=%s%n", rg.getId(), rg.getPlatform(), rg.getLibrary(), rg.getSample()));

    // create the new header and output file
    final SAMFileHeader inHeader = in.getFileHeader();
    final SAMFileHeader outHeader = inHeader.clone();
    outHeader.setReadGroups(Collections.singletonList(rg));
    if (SORT_ORDER != null) outHeader.setSortOrder(SORT_ORDER);

    final SAMFileWriter outWriter = new SAMFileWriterFactory().makeSAMOrBAMWriter(outHeader,
            outHeader.getSortOrder() == inHeader.getSortOrder(),
            OUTPUT);

    final ProgressLogger progress = new ProgressLogger(log);
    for (final SAMRecord read : in) {
        read.setAttribute(SAMTag.RG.name(), RGID);
        outWriter.addAlignment(read);
        progress.record(read);
    }

    // cleanup
    CloserUtil.close(in);
    outWriter.close();
    return 0;
}
 

@Test
public void test1Sample2Readgroups() throws Exception {
    final SAMReadGroupRecord readGroupOne = new SAMReadGroupRecord("rg1");
    readGroupOne.setSample("sample1");
    final SAMReadGroupRecord readGroupTwo = new SAMReadGroupRecord("rg2");
    readGroupTwo.setSample("sample1");

    final SAMFileHeader header = ArtificialReadUtils.createArtificialSamHeader(1, 1, 1000);
    header.addReadGroup(readGroupOne);
    header.addReadGroup(readGroupTwo);

    final Locatable loc = new SimpleInterval("chr1", 1, 1);
    final GATKRead read1 = ArtificialReadUtils.createArtificialRead(header,"read1",0,1,10);
    read1.setReadGroup(readGroupOne.getId());
    final GATKRead read2 = ArtificialReadUtils.createArtificialRead(header,"read2",0,1,10);
    read2.setReadGroup(readGroupTwo.getId());
    final GATKRead read3 = ArtificialReadUtils.createArtificialRead(header,"read3",0,1,10);
    read3.setReadGroup(readGroupOne.getId());
    final GATKRead read4 = ArtificialReadUtils.createArtificialRead(header,"read4",0,1,10);
    read4.setReadGroup(readGroupTwo.getId());
    final GATKRead read5 = ArtificialReadUtils.createArtificialRead(header,"read5",0,1,10);
    read5.setReadGroup(readGroupTwo.getId());
    final GATKRead read6 = ArtificialReadUtils.createArtificialRead(header,"read6",0,1,10);
    read6.setReadGroup(readGroupOne.getId());
    final GATKRead read7 = ArtificialReadUtils.createArtificialRead(header,"read7",0,1,10);
    read7.setReadGroup(readGroupOne.getId());

    final ReadPileup pileup = new ReadPileup(loc, Arrays.asList(read1, read2, read3, read4, read5, read6, read7), 1);

    final AlignmentContext ac = new AlignmentContext(loc, pileup);
    Assert.assertSame(ac.getBasePileup(), pileup);
    Assert.assertEquals(ac.getContig(), loc.getContig());
    Assert.assertEquals(ac.getEnd(), loc.getEnd());
    Assert.assertEquals(ac.getLocation(), loc);
    Assert.assertEquals(ac.getPosition(), loc.getStart());
    Assert.assertEquals(ac.getStart(), loc.getStart());
    Assert.assertEquals(ac.hasPileupBeenDownsampled(), false);
    Assert.assertEquals(ac.size(), pileup.size());

    Assert.assertSame(ac.stratify(AlignmentContext.ReadOrientation.COMPLETE), ac, "Complete");

    final AlignmentContext acFWD = ac.stratify(AlignmentContext.ReadOrientation.FORWARD);
    Assert.assertEquals(acFWD.getLocation(), loc, "Forward Loc");
    Assert.assertEquals((Iterable<?>) acFWD.getBasePileup(), (Iterable<?>)pileup, "Forward Pileup");

    final AlignmentContext acREV = ac.stratify(AlignmentContext.ReadOrientation.REVERSE);
    AlignmentContext emptyAC= new AlignmentContext(loc, new ReadPileup(loc));
    Assert.assertEquals(acREV.getLocation(), loc, "Reverse Loc");
    Assert.assertEquals((Iterable<?>)acREV.getBasePileup(), (Iterable<?>)emptyAC.getBasePileup(), "Reverse pileup");
    Assert.assertNotNull(ac.toString());

    final Map<String, AlignmentContext> bySample = ac.splitContextBySampleName(header);
    Assert.assertEquals(bySample.size(), 1);
    Assert.assertEquals(bySample.keySet(), ReadUtils.getSamplesFromHeader(header));
    final AlignmentContext firstAC = bySample.values().iterator().next();
    Assert.assertEquals(firstAC.getLocation(), ac.getLocation());
    Assert.assertEquals(firstAC.getBasePileup(), ac.getBasePileup());

    final Map<String, AlignmentContext> bySampleAssume1 = ac.splitContextBySampleName("sample1", header);
    Assert.assertEquals(bySampleAssume1.size(), 1);
    Assert.assertEquals(bySampleAssume1.keySet(), ReadUtils.getSamplesFromHeader(header));
    final AlignmentContext firstACAssume1 = bySampleAssume1.values().iterator().next();
    Assert.assertEquals(firstACAssume1.getLocation(), ac.getLocation());
    Assert.assertEquals(firstACAssume1.getBasePileup(), ac.getBasePileup());

    final Map<String, AlignmentContext> stringAlignmentContextMap = AlignmentContext.splitContextBySampleName(pileup, header);
    Assert.assertEquals(stringAlignmentContextMap.keySet(), Collections.singleton("sample1"));
    Assert.assertEquals(stringAlignmentContextMap.get("sample1").getLocation(), loc);
    Assert.assertEquals(stringAlignmentContextMap.get("sample1").getBasePileup(), pileup);
}
 

@Test
public void testOverlappingFragmentFilter() throws Exception {
    final String rgID = "MY.ID";
    final SAMReadGroupRecord rg = new SAMReadGroupRecord(rgID);
    rg.setSample("MySample");
    final SAMFileHeader header = ArtificialReadUtils.createArtificialSamHeaderWithReadGroup(rg);
    final int readlength = 10;
    final String cigar = "10M";
    final byte[] lowQuals = Utils.repeatBytes((byte)10, readlength);
    final byte[] highQuals = Utils.repeatBytes((byte)20, readlength);
    final byte[] basesAllA = Utils.repeatChars('A', readlength);
    final byte[] basesAllC = Utils.repeatChars('C', readlength);
    final byte[] basesAllT = Utils.repeatChars('T', readlength);

    GATKRead read1BasesDisagree = ArtificialReadUtils.createArtificialRead(header, basesAllA, lowQuals, cigar);
    GATKRead read2BasesDisagree = ArtificialReadUtils.createArtificialRead(header, basesAllC, highQuals, cigar);
    GATKRead read3BasesDisagree = ArtificialReadUtils.createArtificialRead(header, basesAllT, lowQuals, cigar);
    read1BasesDisagree.setName("BasesDisagree");
    read2BasesDisagree.setName("BasesDisagree");
    read3BasesDisagree.setName("BasesDisagree");
    read1BasesDisagree.setReadGroup(rgID);
    read2BasesDisagree.setReadGroup(rgID);
    read3BasesDisagree.setReadGroup(rgID);

    GATKRead read1BasesAgree = ArtificialReadUtils.createArtificialRead(header, basesAllA, lowQuals, cigar);
    GATKRead read2BasesAgree = ArtificialReadUtils.createArtificialRead(header, basesAllA, highQuals, cigar);
    read1BasesAgree.setName("BasesAgree");
    read2BasesAgree.setName("BasesAgree");
    read1BasesAgree.setReadGroup(rgID);
    read2BasesAgree.setReadGroup(rgID);

    final List<GATKRead> reads = Arrays.asList(read1BasesDisagree, read2BasesDisagree, read1BasesAgree, read2BasesAgree);
    final int off = 0;
    final ReadPileup pu = new ReadPileup(loc, reads, off);

    final ReadPileup filteredPileupDiscardDiscordant = pu.getOverlappingFragmentFilteredPileup(header);
    final List<GATKRead> filteredReadsDiscardDiscordant = filteredPileupDiscardDiscordant.getReads();
    Assert.assertFalse(filteredReadsDiscardDiscordant.contains(read1BasesDisagree), "Reads with disagreeing bases were kept.");
    Assert.assertFalse(filteredReadsDiscardDiscordant.contains(read2BasesDisagree), "Reads with disagreeing bases were kept.");
    Assert.assertFalse(filteredReadsDiscardDiscordant.contains(read3BasesDisagree), "Reads with disagreeing bases were kept.");
    Assert.assertFalse(filteredReadsDiscardDiscordant.contains(read1BasesAgree), "The lower quality base was kept.");
    Assert.assertTrue(filteredReadsDiscardDiscordant.contains(read2BasesAgree), "The higher quality base is missing.");

    final ReadPileup filteredPileupKeepDiscordant = pu.getOverlappingFragmentFilteredPileup(false, ReadPileup.baseQualTieBreaker, header);
    final List<GATKRead> filteredReadsKeepDiscordant = filteredPileupKeepDiscordant.getReads();
    Assert.assertFalse(filteredReadsKeepDiscordant.contains(read1BasesDisagree), "Low quality base was kept.");
    Assert.assertTrue(filteredReadsKeepDiscordant.contains(read2BasesDisagree), "High quality base is missing.");
    Assert.assertFalse(filteredReadsKeepDiscordant.contains(read3BasesDisagree), "Low quality base was kept.");
    Assert.assertFalse(filteredReadsKeepDiscordant.contains(read1BasesAgree), "The lower quality base was kept.");
    Assert.assertTrue(filteredReadsKeepDiscordant.contains(read2BasesAgree), "The higher quality base is missing.");
}