Unconfigured Ad

Collapse
X
 
  • Time
  • Show
Clear All
new posts
  • Lien
    Member
    • Dec 2009
    • 47

    #1

    Question about python script generating synthetic sequences

    Hi all,

    I'm trying to replicate a Nature Protocols paper: genome-wide copy number analysis of single cells from Baslan et al., 2012.
    In step 69 they generate synthetic 50 bp reads for each position on the hg19 genome: start at position one in the chromosome sequence and take the first 50 bases. Also create read ID strings and quality score strings. These can be Illumina or fastq format. Output these to a file and continue likewise at positions 2, 3, 4 and so on until the end of the chromosome is reached.
    The script to do this can be found in the supplementary data, and I've pasted it below:

    #!/usr/bin/env python

    import re
    import sys
    import random
    import string

    def main():

    CHROMLIST = open("/filepath/chromlist.txt", "r")

    chromlist = []
    for x in CHROMLIST:
    if x.rstrip().find("_") == -1:
    chromlist.append(x.rstrip())

    INFILE = False
    OUTFILE = open("/filepath/sequence.part.0.k50.txt", "w")
    outfilecount = 0
    counter = 0

    for chrFile in chromlist:
    a = chrFile.split("/")
    b = a[len(a) - 1]
    thisChrom = b.split(".")[0]

    if INFILE:
    INFILE.close()
    INFILE = open(chrFile, "r")
    chr = []
    x = ""
    y = ""
    INFILE.readline()
    for x in INFILE:
    chr.append(x.rstrip())
    x = "".join(chr)
    y = x.upper()
    print "after read " + thisChrom
    sys.stdout.flush()

    for i in range(len(y) - 50):
    counter += 1
    thisRead = y[i: (i+50)]

    if counter % 150000000 == 0:
    OUTFILE.close()
    outfilecount += 1
    OUTFILE = open("/filepath/sequence.part." + str(outfilecount) + ".k50.txt", "w")

    OUTFILE.write("@" + thisChrom + "." + str(i))
    OUTFILE.write("\n")
    OUTFILE.write(thisRead)
    OUTFILE.write("\n")
    OUTFILE.write("+" + thisChrom + "." + str(i))
    OUTFILE.write("\n")
    OUTFILE.write("b" * len(thisRead))
    OUTFILE.write("\n")

    INFILE.close()

    OUTFILE.close()


    if __name__ == "__main__":
    main()


    I don't have a lot of experience with python, and when I run this script (after adapting the filepath), I only get 1 empty file. Is there something else I need to adapt?

    Many thanks in advance,
    Lien
    Last edited by GenoMax; 04-18-2016, 02:57 AM. Reason: Moved post to a new thread to increase visibility
  • GenoMax
    Senior Member
    • Feb 2008
    • 7142

    #2
    For reference cross-posted: https://www.biostars.org/p/187267/

    Comment

    Latest Articles

    Collapse

    • SEQadmin2
      Beyond CRISPR/Cas9: Understand, Choose, and Use the Right Genome Editing Tool
      by SEQadmin2



      CRISPR/Cas9 sparked the gene editing revolution for both research and therapeutics.1 But this system still showed severe issues that limited its applications. The most prominent were the heavy reliance on PAM sequences, delivery limitations, double-stranded breaks that prompt unintended edits and cell death, and editing inefficiency (both in targeting and in knock-in reliability).

      Despite this, “CRISPR helped turn genome editing from a specialized technique into
      ...
      07-31-2026, 11:01 AM
    • SEQadmin2
      Proteomic Platforms: How to Choose the Right Analytical Strategy to Improve Detection and Clinical Applications
      by SEQadmin2


      Proteomics platforms are evolving rapidly, with advances in mass spectrometry and affinity-based approaches expanding what researchers can detect and at what scale. As the field moves toward deeper proteome coverage and clinical applications, scientists face an increasingly complex landscape of tools. This article will explore how researchers are navigating these choices to find the right platform for their work.

      The systematic characterization of the human proteome has
      ...
      07-20-2026, 11:48 AM

    ad_right_rmr

    Collapse

    News

    Collapse

    Topics Statistics Last Post
    Started by SEQadmin2, 08-13-2026, 12:22 PM
    0 responses
    31 views
    0 reactions
    Last Post SEQadmin2  
    Started by SEQadmin2, 08-11-2026, 10:35 AM
    0 responses
    24 views
    0 reactions
    Last Post SEQadmin2  
    Started by SEQadmin2, 08-06-2026, 07:41 AM
    0 responses
    38 views
    0 reactions
    Last Post SEQadmin2  
    Started by SEQadmin2, 08-03-2026, 10:13 AM
    0 responses
    51 views
    0 reactions
    Last Post SEQadmin2  
    Working...