使用telomerehunter2计算端粒长度

官方github链接：https://github.com/ferdinand-popp/telomerehunter2

参考代码：

telomerehunter2 -ibt TUMOR_FILE -ibc CONTROL_FILE -o OUTPUT_DIRECTORY -p ID_OF_SAMPLE -b BANDING_FILE [options]

批量运行脚本：

perl -ne 'chomp; next if /^$/; @a = split /\t/; print "telomerehunter2 -p $a[0] -ibc ../align/$a[1]_bqsr.bam -ibt ../align/$a[0]_bqsr.bam -o /group_homes/prolinying/home/share/Phyllodes/LCM-WGS/telomerehunter2 -b /home/zhoukaiwen/database/telomerehunter2/hg38_cytoBand.txt > log/$a[0].log 2>&1 && echo $a[0] telomerehunter2 ok\n";' ../mutect2/sample_pair.txt > RunTelomereHunter.sh

nohup bash -c "cat RunTelomereHunter.sh | parallel -j 4" > log/telomerehunter.log 2>&1 &

结果文件解读：
输出文件为“Tumor_summary.tsv”的即为结果，其中tel_content为端粒长度

PID     sample  tel_content     total_reads     read_lengths    repeat_threshold_set    repeat_threshold_used   intratelomeric_reads    junctionspanning_reads  subtelomeric_reads      intrachromosomal_reads  tel_read_countgc_bins_for_correction  total_reads_with_tel_gc TCAGGG_arbitrary_context_norm_by_intratel_reads TGAGGG_arbitrary_context_norm_by_intratel_reads TTGGGG_arbitrary_context_norm_by_intratel_reads TTCGGG_arbitrary_context_norm_by_intratel_reads       TTTGGG_arbitrary_context_norm_by_intratel_reads ATAGGG_arbitrary_context_norm_by_intratel_reads CATGGG_arbitrary_context_norm_by_intratel_reads CTAGGG_arbitrary_context_norm_by_intratel_reads GTAGGG_arbitrary_context_norm_by_intratel_reads       TAAGGG_arbitrary_context_norm_by_intratel_reads TCAGGG_singletons_norm_by_all_reads     TGAGGG_singletons_norm_by_all_reads     TTGGGG_singletons_norm_by_all_reads     TTCGGG_singletons_norm_by_all_reads   TTTGGG_singletons_norm_by_all_reads     ATAGGG_singletons_norm_by_all_reads     CATGGG_singletons_norm_by_all_reads     CTAGGG_singletons_norm_by_all_reads     GTAGGG_singletons_norm_by_all_reads   TAAGGG_singletons_norm_by_all_reads
FETB06-BLPT-M   tumor   2704.045416     2376228562      90      6 per 100 bp    5       1116591 13770   144649  53817   1328827 48-52   412933523       0.0943326607504448      0.1466427725102566      0.083708358745503    0.0422858504143415       0.0132438824959183      0.0180191314456233      0.0006931812991507      0.0899380346071211      0.0407203712012724      0.0246464461920255      4.965852270569585e-07   7.01111007014316e-07    4.717559657041106e-07 2.2725086661928608e-08  7.221527539235092e-07   3.577096974562837e-07   0       5.260436727298289e-07   1.658089656444421e-07   1.392963645388587e-07
FETB06-BLPT-M   control 1775.935218     2160559050      90      6 per 100 bp    5       491170  6264    47324   23042   567800  48-52   276569773       0.0725003562921188      0.1134841297310503      0.0863672455565282   0.0403058004356943       0.0168577071075187      0.0165054868986298      0.0011767819695828      0.0497526314717918      0.0713927967913349      0.0369810859783781      3.3833835738023454e-07  2.615063911352018e-07   3.813827722042589e-07 3.7027453612063976e-08  4.670087586821568e-07   8.192324111669153e-08   0       5.688342561153328e-07   1.383901078750891e-07   2.2077619216193147e-07
FETB06-BLPT-M   log2(tumor/control)     0.6065404262406129                                                                                              0.37976927592854426     0.3698154107408493      -0.0451125885059767140.06918752014628017      -0.34808220635008424    0.12658377480988212     -0.7635424085069049     0.8541585240923169      -0.8100278073435476     -0.5854079539855125     0.5535745987608439      1.422796678527075       0.306801086086286     -0.7043095422598166     0.6288544288572661      2.126444531406144               -0.11282577035913463    0.2607811957223489      -0.6644269987342165

合并输出结果：
merge_summary.py

import pandas as pd
import os
import glob

def merge_telomerehunter_outputs():
    # 1. 获取当前目录下所有样本文件夹中的 summary.tsv 文件路径
    # 假设文件夹结构为: 样本名/样本名_summary.tsv
    pattern = "*/*_summary.tsv"
    files = glob.glob(pattern)

    if not files:
        print("未找到任何 *_summary.tsv 文件，请检查当前目录结构。")
        return

    all_data = []

    # 2. 遍历所有文件并读取
    for file in files:
        try:
            # telomerehunter2 输出通常为制表符分隔，这里使用 sep='\t'
            df = pd.read_csv(file, sep='\t')
            all_data.append(df)
        except Exception as e:
            print(f"读取文件 {file} 时出错: {e}")

    if not all_data:
        print("没有成功读取到任何有效数据。")
        return

    # 3. 合并所有数据
    combined_df = pd.concat(all_data, ignore_index=True)

    # 4. 根据 'sample' 列拆分为三个文件
    # 注意：这里假设样本类型严格为 'tumor', 'control', 'log2(tumor/control)'
    tumor_df = combined_df[combined_df['sample'] == 'tumor']
    normal_df = combined_df[combined_df['sample'] == 'control']
    log2_df = combined_df[combined_df['sample'] == 'log2(tumor/control)']

    # 5. 导出为三个独立的 TSV 文件
    tumor_df.to_csv('merged_tumor_summary.tsv', sep='\t', index=False)
    normal_df.to_csv('merged_normal_summary.tsv', sep='\t', index=False)
    log2_df.to_csv('merged_log2_summary.tsv', sep='\t', index=False)

    print(f"处理完成！共整合了 {len(files)} 个样本。")
    print(f"- Tumor 样本数: {len(tumor_df)}")
    print(f"- Normal 样本数: {len(normal_df)}")
    print(f"- Log2 样本数: {len(log2_df)}")

if __name__ == "__main__":
    merge_telomerehunter_outputs()