Commit 16bc0e80 authored by Jakob Zierk's avatar Jakob Zierk

Initial commit.

parent 33964261
.vscode/
.env
.pyc
anonymous-patient-data/
from kosmic import kosmic
import json
rs = {}
for analyte in ["Hb", "Leu", "Thr"]:
for group in ["A", "AB", "ABC", "ABCD", "BC", "BCD", "CD"]:
decimals = { "Hb": 1, "Leu": 1, "Thr": 0 }[analyte]
with open("anonymous-patient-data/" + analyte + "-" + group + ".txt") as f:
samples = [float(line) for line in f]
rs[analyte + "-" + group] = kosmic(samples, decimals, bootstrap=100)
with open("results/patient-data.json", "w") as json_file:
json.dump(rs, json_file)
from kosmic import kosmic, percentile
import random
import numpy as np
from multiprocessing import Pool, cpu_count
# Simulation:
mu0 = 14.0
n = 10000
def sim(mu1, mu2, pi1, pi2):
random.seed(0)
p025s = []
p975s = []
for _ in range(100):
samples = []
for _ in range(int(n * (1-pi1-pi2))):
samples.append(round(random.normalvariate(mu0, 1/0.979982), 1))
for _ in range(int(n * pi1)):
samples.append(round(random.normalvariate(mu1, 1/0.979982), 1))
for _ in range(int(n * pi2)):
samples.append(round(random.normalvariate(mu2, 1/0.979982), 1))
random.shuffle(samples)
result = kosmic(samples, 1)
p025s.append(percentile(result, 0.025))
p975s.append(percentile(result, 0.975))
return p025s, p975s
if __name__ == "__main__":
results = {}
with Pool(processes=cpu_count()) as pool:
# Start simulations:
for mu1 in [15, 16, 17, 18]:
for mu2 in [10, 11, 12, 13]:
for pi1 in [0.0, 0.05, 0.1, 0.15, 0.2, 0.3]:
for pi2 in [0.0, 0.05, 0.1, 0.15, 0.2, 0.3]:
if pi1 == 0.0 and pi2 == 0.0 and (mu1 != 15.0 or mu2 != 13.0):
continue
results[str(mu1) + "\t" + str(mu2) + "\t" + str(pi1) + "\t" + str(pi2)] = pool.apply_async(sim, [mu1, mu2, pi1, pi2])
# Output simulation results:
with open("results/hemoglobin-simulation.tsv", "w", encoding="utf-8") as hb:
hb.write("\tµ₁\tµ₂\tπ₁\tπ₂\tPᵀ\t\tT\t\t\t\n")
for task, result in results.items():
result = result.get()
hb.write("Median\t" + task + "\t" + str(np.percentile(result[0], 50)) + "\t" + str(np.percentile(result[1], 50)) + "\n")
hb.write("P5\t" + task + "\t" + str(np.percentile(result[0], 5)) + "\t" + str(np.percentile(result[1], 5)) + "\n")
hb.write("P95\t" + task + "\t" + str(np.percentile(result[0], 95)) + "\t" + str(np.percentile(result[1], 95)) + "\n")
import random
import numpy as np
from multiprocessing import Pool, cpu_count
from kosmic import kosmic, percentile
# Simulation:
mu0 = 14.0
n = 10000
def sim(mu1, mu2, pi1, pi2):
random.seed(0)
p025s = []
p975s = []
for _ in range(100):
samples = []
for _ in range(int(n * (1-pi1-pi2))):
samples.append(round(random.normalvariate(mu0, 1/0.979982), 1))
for _ in range(int(n * pi1)):
samples.append(round(random.normalvariate(mu1, 1/0.979982), 1))
for _ in range(int(n * pi2)):
samples.append(round(random.normalvariate(mu2, 1/0.979982), 1))
random.shuffle(samples)
result = kosmic(samples, 1)
p025s.append(percentile(result, 0.025))
p975s.append(percentile(result, 0.975))
return p025s, p975s
with open("results/publication-figure-2-b.tsv", "w", encoding="utf-8") as fig2:
for pi1 in np.arange(0.0, 0.51, 0.01):
result = sim(10, 18, pi1, 0)
fig2.write(str(pi1) + "\t" + str(np.percentile(result[0], 50)) + "\t" + str(np.percentile(result[0], 5)) + "\t" + str(np.percentile(result[0], 95)) + "\n")
with open("results/publication-figure-2-c.tsv", "w", encoding="utf-8") as fig2:
for mu1 in np.arange(8.0, 14.1, 0.1):
result = sim(mu1, 18, 0.20, 0)
fig2.write(str(mu1) + "\t" + str(np.percentile(result[0], 50)) + "\t" + str(np.percentile(result[0], 5)) + "\t" + str(np.percentile(result[0], 95)) + "\n")
from kosmic import kosmic, percentile
import random
import numpy as np
from multiprocessing import Pool, cpu_count
from math import exp
meanlog=-5.506435e-06
sdlog=7.073083e-01
# Simulation:
n = 50000
decimals = 2
def sim(mu1, sd1, pi1):
random.seed(0)
samples = []
for _ in range(int(n * (1-pi1))):
samples.append(round(random.lognormvariate(-5.506435e-06, 7.073083e-01), decimals))
for _ in range(int(n * pi1)):
samples.append(round(random.normalvariate(mu1, sd1/0.979982), decimals))
random.shuffle(samples)
result = kosmic(samples, decimals)
return percentile(result, 0.025), percentile(result, 0.975)
if __name__ == "__main__":
results = {}
with Pool(processes=cpu_count()) as pool:
# Start simulations:
for mu1 in [3, 3.5, 4, 4.5, 5, 5.5, 6]:
for sd1 in [.5, 1, 1.5, 2]:
for pi1 in [0.0, 0.05, 0.1, 0.15, 0.2, 0.3]:
if pi1 == 0.0 and (mu1 != 3 or sd1 != .5):
continue
results[str(mu1) + "\t" + str(sd1) + "\t" + str(pi1)] = pool.apply_async(sim, [mu1, sd1, pi1])
# Output simulation results:
with open("results/tsh-simulation.tsv", "w", encoding="utf-8") as tsh:
tsh.write("\tµ₁\tsd₁\tπ₁\tPᵀ\t\n")
for task, result in results.items():
result = result.get()
tsh.write("Result\t" + task + "\t" + str(result[0]) + "\t" + str(result[1]) + "\n")
from kosmic import kosmic, percentile
import random
import numpy as np
from multiprocessing import Pool, cpu_count
from math import exp
# Simulation:
n = 25000
decimals = 0
def sim(a1, pi1):
random.seed(0)
p025s = []
p975s = []
for _ in range(100):
samples = []
for _ in range(int(n * (1-pi1))):
samples.append(round(random.lognormvariate(3.1073040, 0.4105784), decimals))
for _ in range(int(n * (pi1))):
samples.append(round(a1*random.lognormvariate(3.1073040, 0.4105784), decimals))
random.shuffle(samples)
result = rle2(samples, decimals)
p025s.append(percentile(result, 0.025))
p975s.append(percentile(result, 0.975))
return p025s, p975s
if __name__ == "__main__":
results = {}
with Pool(processes=cpu_count()) as pool:
# Start simulations:
for a1 in [.1, .25, .5, 1.5, 2.0, 2.5, 3.0, 5.0]:
for pi1 in [0.0, 0.05, 0.1, 0.15, 0.2, 0.3]:
if pi1 == 0.0 and (a1 != .1):
continue
results[str(a1) + "\t" + str(pi1)] = pool.apply_async(sim, [a1, pi1])
# Output simulation results:
with open("results/y-gt-simulation.tsv", "w", encoding="utf-8") as ygt:
ygt.write("\ta₁\tπ₁\tPᵀ\t\n")
for task, result in results.items():
result = result.get()
ygt.write("Median\t" + task + "\t" + str(np.percentile(result[0], 50)) + "\t" + str(np.percentile(result[1], 50)) + "\n")
ygt.write("P5\t" + task + "\t" + str(np.percentile(result[0], 5)) + "\t" + str(np.percentile(result[1], 5)) + "\n")
ygt.write("P95\t" + task + "\t" + str(np.percentile(result[0], 95)) + "\t" + str(np.percentile(result[1], 95)) + "\n")
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment