# считаем метрики TP/FP/FN/TN для разных порогов hmmsearch
# true positives - reviewed-белки с обоими доменами PF00072 + PF03861

tp_ids = set(open('true_positives.txt').read().split())
total = sum(1 for line in open('family_reviewed.fasta') if line.startswith('>'))

print(f"всего белков в области поиска: {total}")
print(f"true positives: {len(tp_ids)}")

# читаем скоры из tableoutput.txt
# имена в формате sp|ID|name, скор в колонке 6 (индекс 5)
scores = {}
with open('tableoutput.txt') as f:
    for line in f:
        if line.startswith('#') or not line.strip():
            continue
        parts = line.split()
        name = parts[0]
        uid = name.split('|')[1] if '|' in name else name
        score = float(parts[5])
        scores[uid] = score

print(f"нашлось хитов: {len(scores)}")
print()

# для каждого уникального значения скора считаем метрики
thresholds = sorted(set(scores.values()), reverse=True)

print(f"{'порог':>8}  {'TP':>4}  {'FP':>6}  {'FN':>4}  {'TN':>6}  {'сумма':>6}")
for t in thresholds:
    TP = sum(1 for uid, s in scores.items() if uid in tp_ids and s >= t)
    FP = sum(1 for uid, s in scores.items() if uid not in tp_ids and s >= t)
    FN = len(tp_ids) - TP
    TN = total - TP - FP - FN
    print(f"{t:>8.1f}  {TP:>4}  {FP:>6}  {FN:>4}  {TN:>6}  {TP+FP+FN+TN:>6}")