Here’s the modeling code from the substack post.<\/p>\n\n\n\n
from math import comb\nfrom math import log10\nimport sys\n\n\ndef tprb(n, k, t):\n\n totalp = 1\n\n for i in range(0,k,1):\n p = 1 - (t \/ (n - i))\n totalp = totalp * p\n\n return 1-totalp\n\nctoffset = 45\nctslope = -3.36\n\ndef predictCT(fragments,ampregion,fragmentsize,genomesize):\n\ttotalfragments = genomesize-fragmentsize\n\tontargetfragments = fragmentsize-ampregion\n\ttargetfragments = fragments*(ontargetfragments\/totalfragments)\n\tct = ctoffset+log10(targetfragments)*ctslope\n\treturn ct\n\ngenomesize = 30000\ntotalng = 10\nqpcrampregion = 90\nbasesperng = 1771000000000 \nrrnafraction = 0.6\nfragmentsize = 1000\nstarttargetfraction = 0.000000001\ntargetfraction = starttargetfraction\n\n\ntotalfragments = (basesperng*totalng)\/fragmentsize\nreads = 100000\ntargetfragments = totalfragments*targetfraction\n\ntotalfragments_norRNA = ((totalfragments-targetfragments)*(1-rrnafraction)) + targetfragments\n\n# Vary reads count\n#print (\"total fragments, reads, target fragments, probability of >1 read, probability of >1 read (no rRNA), predicted CT\")\n#for creads in range(reads,100000000,1000000):\n#\tprint (totalfragments, \" \", creads, \" \", targetfragments, \" \", tprb(totalfragments,creads,targetfragments), \" \", tprb(totalfragments_norRNA,creads,targetfragments),predictCT(targetfragments,qpcrampregion,fragmentsize,genomesize))\n\n# Vary target fraction\n#print (\"target fraction, total fragments, reads, target fragments, probability of >1 read, probability of >1 read (no rRNA), predicted CT\")\n#while targetfraction < 1:\n#\ttargetfragments = totalfragments*targetfraction\n#\tprint (targetfraction, \" \", totalfragments, \" \", reads, \" \", targetfragments, \" \", tprb(totalfragments,reads,targetfragments), \" \", tprb(totalfragments_norRNA,reads,targetfragments),predictCT(targetfragments,qpcrampregion,fragmentsize,genomesize))\n#\ttargetfraction *= 10\n\n# Vary both\n#print(\", \", end='')\n#while targetfraction < 1:\n#\ttargetfragments = totalfragments*targetfraction\n#\tprint(targetfraction, \", \", end='')\n#\ttargetfraction *= 2\n#print()\n\n#for creads in range(reads,5000000,100000):\n#\tprint(creads, \", \", end='')\n#\ttargetfraction = starttargetfraction\n#\twhile targetfraction < 1:\n#\t\ttargetfragments = totalfragments*targetfraction\n#\t\tprint (tprb(totalfragments,creads,targetfragments), \", \", end='')\n#\t\t\n#\t\ttargetfraction *= 2\n#\tprint()\n\n# Vary both - CT\nprint(\", \", end='')\nwhile targetfraction < 1:\n\ttargetfragments = totalfragments*targetfraction\n\tprint(\"{:.2f}\".format(predictCT(targetfragments,qpcrampregion,fragmentsize,genomesize)), \", \", end='')\n\ttargetfraction *= 2\nprint()\n\nfor creads in range(reads,5000000,100000):\n\tprint(creads, \", \", end='')\n\ttargetfraction = starttargetfraction\n\twhile targetfraction < 1:\n\t\ttargetfragments = totalfragments*targetfraction\n\t\tprint (tprb(totalfragments,creads,targetfragments), \", \", end='')\n\t\t\n\t\ttargetfraction *= 2\n\tprint()\n<\/code><\/pre>\n","protected":false},"excerpt":{"rendered":"Here’s the modeling code from the substack post.<\/p>\n","protected":false},"author":1,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":{"jetpack_post_was_ever_published":false,"_jetpack_newsletter_access":"","_jetpack_dont_email_post_to_subs":false,"_jetpack_newsletter_tier_id":0,"_jetpack_memberships_contains_paywalled_content":false,"_jetpack_memberships_contains_paid_content":false,"footnotes":"","jetpack_publicize_message":"","jetpack_publicize_feature_enabled":true,"jetpack_social_post_already_shared":true,"jetpack_social_options":{"image_generator_settings":{"template":"highway","default_image_id":0,"font":"","enabled":false},"version":2}},"categories":[1],"tags":[],"class_list":["post-6788","post","type-post","status-publish","format-standard","hentry","category-uncategorized"],"jetpack_publicize_connections":[],"jetpack_featured_media_url":"","jetpack_shortlink":"https:\/\/wp.me\/p1RRoU-1Lu","jetpack_sharing_enabled":true,"_links":{"self":[{"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/posts\/6788","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/users\/1"}],"replies":[{"embeddable":true,"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/comments?post=6788"}],"version-history":[{"count":1,"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/posts\/6788\/revisions"}],"predecessor-version":[{"id":6789,"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/posts\/6788\/revisions\/6789"}],"wp:attachment":[{"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/media?parent=6788"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/categories?post=6788"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/41j.com\/blog\/wp-json\/wp\/v2\/tags?post=6788"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}