根据基因列表，统计SV注释文件中，基因受到SV影响的数目

import sys
def read_gene_list(filename):
    """从文件中读取基因列表"""
    with open(filename, 'r') as file:
        genes = [line.strip() for line in file]
    return genes
def process_sv_file(genes, sv_filename, output_filename):
    """处理 SV 文件并筛选出基因列表中的基因及其对应的 SV 类型"""
    sv_type_order = ['INS', 'DEL', 'DUP', 'INV']
    gene_sv_counts = {gene: {sv_type: 0 for sv_type in sv_type_order} for gene in genes}
    with open(sv_filename, 'r') as sv_file:
        for line in sv_file:
            parts = line.strip().split('\t')
            if len(parts) >= 11:
                gene_part = parts[6]
                for gene in genes:
                    if gene in gene_part:
                        sv_type = parts[10]
                        if sv_type in sv_type_order:
                            gene_sv_counts[gene][sv_type] += 1
                        break

    with open(output_filename, 'w') as output_file:
        # 按照基因列表的顺序输出所有基因及其每种 SV 类型的计数
        for gene in genes:
            output_file.write(f"{gene}\t")
            for sv_type in sv_type_order:
                output_file.write(f"{sv_type}:{gene_sv_counts[gene][sv_type]}\t")
            output_file.write("\n")


if __name__ == "__main__":
    if len(sys.argv)!= 4:
        print("Usage: python script.py <gene_list_file> <sv_data_file> <output_file>")
        sys.exit(1)

    gene_list_filename = sys.argv[1]
    sv_data_filename = sys.argv[2]
    output_filename = sys.argv[3]

    # 读取基因列表
    genes = read_gene_list(gene_list_filename)

    # 处理 SV 数据文件
    process_sv_file(genes, sv_data_filename, output_filename)

    print("处理完成！")

cat all|while read id
do
python wo4.py gene ../$id add1/$id
done

结果

合并表格

import os
import csv

def merge_tsv_files(output_filename):
    """合并当前目录下的所有 TSV 文件到一个新的 TSV 文件中，并保留文件名作为额外一列"""
    
    # 获取当前目录下的所有文件
    files = os.listdir('.')
    
    # 过滤出 TSV 文件
    tsv_files = [f for f in files if f.endswith('.tsv')]
    
    # 如果没有找到任何 TSV 文件，则返回
    if not tsv_files:
        print("没有找到 TSV 文件。")
        return
    
    # 写入头部信息（假设所有文件的头部都相同）
    header_written = False

    with open(output_filename, 'w', newline='') as outfile:
        writer = csv.writer(outfile, delimiter='\t')

        for filename in tsv_files:
            with open(filename, 'r') as infile:
                reader = csv.reader(infile, delimiter='\t')
                
                # 只写入一次头部信息
                if not header_written:
                    header = next(reader)
                    # 添加文件名列
       #             header.append('File_Name')
        #            writer.writerow(header)
                    header_written = True
                
                # 写入数据行
                for row in reader:
                    # 添加文件名作为最后一列
                    row.append(filename)
                    writer.writerow(row)
    
    print(f"文件已成功合并到 {output_filename}")

if __name__ == "__main__":
    # 指定输出文件名
    output_filename = "merged_with_filenames.tsv"
    
    # 合并文件
    merge_tsv_files(output_filename)

R统计结果

setwd("D:\\科研助理工作\\data\\10-11")

# 读取第一个表格
group_table <- read.table("group.tsv", header = TRUE, sep = "\t")

# 读取第二个表格
merged_table <- read.table("merged_with_filenames.tsv", header = TRUE, sep = "\t")

# 按照'sample'列进行合并
combined_table <- merge(group_table, merged_table, by = "sample")

# 将合并后的表格保存为新的文件
write.table(combined_table, "combined_output.tsv", sep = "\t", row.names = FALSE)



# 读取文件
data <- read.table("combined_output.tsv", header = TRUE, sep = "\t")

# 安装并加载 plyr 包（如果未安装）
if (!require("plyr")) {
  install.packages("plyr")
  library(plyr)
}

# 按照 k5 和 gene 列分组计算均值和标准差，结果保留两位小数
result <- ddply(data,.(k5, gene), summarise,
                INS_mean = round(mean(INS, na.rm = TRUE), 2),
                INS_sd = round(sd(INS, na.rm = TRUE), 2),
                DEL_mean = round(mean(DEL, na.rm = TRUE), 2),
                DEL_sd = round(sd(DEL, na.rm = TRUE), 2),
                DUP_mean = round(mean(DUP, na.rm = TRUE), 2),
                DUP_sd = round(sd(DUP, na.rm = TRUE), 2),
                INV_mean = round(mean(INV, na.rm = TRUE), 2),
                INV_sd = round(sd(INV, na.rm = TRUE), 2),
                total_mean = round(mean(total, na.rm = TRUE), 2),
                total_sd = round(sd(total, na.rm = TRUE), 2))

# 将结果保存为新的表格
write.table(result, "result_table.tsv", sep = "\t", row.names = FALSE)




# 读取文件
data <- read.table("combined_output.tsv", header = TRUE, sep = "\t")

# 安装并加载 plyr 包（如果未安装）
if (!require("plyr")) {
  install.packages("plyr")
  library(plyr)
}

# 按照 k5 和 gene 列分组计算均值和标准差
result <- ddply(data,.(k5, gene), summarise,
                INS_mean = mean(INS, na.rm = TRUE),
                INS_sd = sd(INS, na.rm = TRUE),
                DEL_mean = mean(DEL, na.rm = TRUE),
                DEL_sd = sd(DEL, na.rm = TRUE),
                DUP_mean = mean(DUP, na.rm = TRUE),
                DUP_sd = sd(DUP, na.rm = TRUE),
                INV_mean = mean(INV, na.rm = TRUE),
                INV_sd = sd(INV, na.rm = TRUE),
                total_mean = mean(total, na.rm = TRUE),
                total_sd = sd(total, na.rm = TRUE))

# 将结果保存为新的表格
write.table(result, "result_table.tsv", sep = "\t", row.names = FALSE)






# 读取文件
data <- read.table("combined_output.tsv", header = TRUE, sep = "\t")

# 安装并加载 plyr 包（如果未安装）
if (!require("plyr")) {
  install.packages("plyr")
  library(plyr)
}

# 定义一个函数来格式化均值和标准差为均值±标准差的形式
format_mean_sd <- function(mean_value, sd_value) {
  return(paste0(round(mean_value, 2), "±", round(sd_value, 2)))
}

# 按照 k5 和 gene 列分组计算均值和标准差，并格式化为均值±标准差
result <- ddply(data,.(k5, gene), summarise,
                INS = format_mean_sd(mean(INS, na.rm = TRUE), sd(INS, na.rm = TRUE)),
                DEL = format_mean_sd(mean(DEL, na.rm = TRUE), sd(DEL, na.rm = TRUE)),
                DUP = format_mean_sd(mean(DUP, na.rm = TRUE), sd(DUP, na.rm = TRUE)),
                INV = format_mean_sd(mean(INV, na.rm = TRUE), sd(INV, na.rm = TRUE)),
                total = format_mean_sd(mean(total, na.rm = TRUE), sd(total, na.rm = TRUE)))

# 将结果保存为新的表格
write.table(result, "result_table3.tsv", sep = "\t", row.names = FALSE, quote = FALSE)

结果