计算每个基因的平均覆盖率答案

【问题标题】：Calculating mean coverage per gene计算每个基因的平均覆盖率
【发布时间】：2021-07-23 00:50:56
【问题描述】：

我有 2 个文件：文件 1（下）有 bp 开始和停止坐标

Chrom	Gene start (bp)	Gene end (bp)
1	50902700	50902978
1	103817769	103828355

文件 2 具有每个碱基对位置的平均覆盖率值：

Chrom	pos	mean
1	12141	0.029005
1	12142	0.029216

我需要什么：我需要匹配 chrom, start 和 end 从文件 1（以 start 和 end 作为范围）到 chrom, pos；并计算文件 1 坐标范围内的平均值（文件 2 中的平均列）。

期望的输出：染色体/支架名称基因起点 (bp) 基因终点 (bp) 每个基因的平均覆盖率

Chrom	Gene start (bp)	Gene end (bp)	average of means
1	50902700	50902978	(mean coverage for this gene)
1	103817769	103828355	(mean coverage for this gene)

我尝试过使用字典和嵌套的 for 循环：

代码：

`# importing gene start/end files
df_gene = pd.read_csv('gene_list.csv')
# importing exome data file
df_data = pd.read_csv('exomes.coverage.summary.tsv', sep = '\t')
# Creating a Dictionary to store mean values
dict_mean=df_data.set_index('pos')['mean'].to_dict()

start = df_gene['Gene start (bp)'].to_list()
end = df_gene['Gene end (bp)'].to_list()

list_mean=[]
x=0
df_mean = pd.DataFrame(columns=['start','end','mean coverage'])
### looping:
for s,e in zip(start,end):    
    for key,val in dict_mean.items():
        if key>=s and key<=e:
                list_mean.append(val)
                x=np.mean(list_mean)    #calculating average of means           
                
    my_series = pd.Series(data=[s, e, x], index=['start', 'end', 'mean coverage'])
    df_mean=df_mean.append(my_series,ignore_index=True)
    
### Add mean coverage to gene dataframe
df_gene['mean coverage'] = df_mean['mean coverage']

df_gene.to_csv('gene_out.csv', index=False)
`

该代码有效，但它不考虑色度数。如何找到 start 和 stop 范围内的平均值？

【问题讨论】：

您好，欢迎来到 SO。恐怕我对生物学/生物信息学的了解不够，无法理解这个问题。您能否更清楚地描述您要执行哪些操作，以及“覆盖”代表什么？
表格不完整，有更多数据点，我正在尝试将“chrom”、“start”和“end”作为一个范围，并从文件 2 中找到平均值（文件 2具有范围内各个点的平均值） P.S：覆盖率代表该位置的“排序”程度。谢谢！
请edit您的问题并添加您尝试过的代码并描述它如何不适合您。见How to Ask。

标签： python pandas dataframe bioinformatics

【解决方案1】：

这应该可以解决您的问题。它非常依赖于具有一致的数据，这意味着文件 2 中的每个条目 (chrom, pos) 在文件 1 中都有对应的 (chrom, s, e)。如果不是这种情况，则必须在内部 while 中执行附加检查循环。

# importing gene start/end files
df_gene = pd.read_csv('gene_list.csv')
# importing exome data file
df_data = pd.read_csv('exomes.coverage.summary.tsv', sep = '\t')
# Creating a Dictionary to store mean values
chroms_f2 = df_data.['Chrom'].to_list()
positions = df_data.['pos'].to_list()
means = df_data.['mean'].to_list()
f2_as_list = sorted(zip(chroms_f2, zip(positions, means))

starts = df_gene['Gene start (bp)'].to_list()
ends = df_gene['Gene end (bp)'].to_list()
chroms_f1 = df_gene['Chrom'].to_list()
f1_as_list = sorted(zip(chroms_f1, zip(starts, ends)))

df_mean = pd.DataFrame(columns=['chrom', 'start','end','mean coverage'])

### looping:
i1 = 0
c1, (s, e) = f1_as_list[i1]
list_mean = []
for c2, (p, m) in f2_as_list:
    if not (c1 == c2 and s <= p <= e):
        my_series = pd.Series(
            data=[c, s, e, np.mean(list_mean)], 
            index=['chrom', 'start', 'end', 'mean coverage']
        )
        df_mean=df_mean.append(my_series,ignore_index=True)
        list_mean = []
    while not (c1 == c2 and s <= p <= e):
        i1 += 1
        c1, (s, e) = f1_as_list[i1]
    list_mean.append(m)
    
### Add mean coverage to gene dataframe
df_gene['mean coverage'] = df_mean['mean coverage']

df_gene.to_csv('gene_out.csv', index=False)

【讨论】：