code: def wgsim(args):
      
      p = OptionParser(wgsim.__doc__)
      p.add_option("--erate", default=.01, type="float",
                   help="Base error rate of the read [default: %default]")
      p.add_option("--noerrors", default=False, action="store_true",
                   help="Simulate reads with no errors [default: %default]")
      p.add_option("--genomesize", type="int",
                   help="Genome size in Mb [default: estimate from data]")
      add_sim_options(p)
      opts, args = p.parse_args(args)
  
      if len(args) != 1:
          sys.exit(not p.print_help())
  
      fastafile, = args
      pf = op.basename(fastafile).split(".")[0]
  
      genomesize = opts.genomesize
      size = genomesize * 1000000 if genomesize else Fasta(fastafile).totalsize
      depth = opts.depth
      readlen = opts.readlen
      readnum = int(math.ceil(size * depth / (2 * readlen)))
  
      distance = opts.distance
      stdev = distance / 10
  
      outpf = opts.outfile or "{0}.{1}bp.{2}x".format(pf, distance, depth)
  
      logging.debug("Total genome size: {0} bp".format(size))
      logging.debug("Target depth: {0}x".format(depth))
      logging.debug("Number of read pairs (2x{0}): {1}".format(readlen, readnum))
  
      if opts.noerrors:
          opts.erate = 0
  
      cmd = "dwgsim -e {0} -E {0}".format(opts.erate)
      if opts.noerrors:
          cmd += " -r 0 -R 0 -X 0 -y 0"
  
      cmd += " -d {0} -s {1}".format(distance, stdev)
      cmd += " -N {0} -1 {1} -2 {1}".format(readnum, readlen)
      cmd += " {0} {1}".format(fastafile, outpf)
      sh(cmd)
docs: %prog wgsim fastafile
  
      Run dwgsim on fastafile.