como dividir valores em cada linha pelos valores de soma de duas linhas (traduzir bash para perl ou python)

Question

Eu tenho uma população de dados de contagem de alelos, parece com isso

1   0   0   0   0   0   0   0   0   0   1   2   1   0   0   0   0
0   2   0   0   0   0   0   0   0   0   0   4   0   2   0   0   0
0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
2   2   0   0   2   1   0   0   0   0   2   4   0   0   0   2   0

Colunas são populações e linhas são SNPs. Existem 2 linhas por SNP (uma linha para o número de cópias do alelo "1" em cada população e uma linha para o número de cópias do alelo "2"). No exemplo acima, a primeira e a segunda linha são o número de alelos 1 e 2 para o SNP1, a terceira e a quarta linha são alelos para o SNP 2 e, da mesma forma, para todo o conjunto de dados. Quero calcular a frequência do alelo da população por SNP para todas as populações: frequency of allele 1 at SNP1 in population 1 =Number of copies of allele 1 in population/Total number of 1+2 alleles copies in population e frequency of allele 2 at SNP1 in population 1 =Number of copies of allele a in population/Total number of 1+2 gene copies in population  Isto significa que para cada número SNP de cópias de cada alelo deve ser dividido a soma do número de cópias do alelo "1" e "2" para cada população. Esta é a saída desejada:

1   0   0   0   0   0   0   0   0   0   1   0.333333    1   0   0   0   0
0   1   0   0   0   0   0   0   0   0   0   0.666667    0   1   0   0   0
0   0   0   0   0   0   0   0   0   0   0   0           0   0   0   0   0
1   1   0   0   1   1   0   0   0   0   1   1           0   0   0   1   0

Eu tenho R e bash solução, mas existe alguma maneira de fazer essa estimativa em perl ou python? Qualquer ajuda é apreciada.

Eu tenho a solução bash, aqui

awk '{for(i=1; i<=NF; i++) tally[i]+=$i}; (NR%2)==1{for(i=1; i<=NF; i++) allele1[i]=$i}; (NR%2)==0{for(i=1; i<=NF; i++) allele2[i]=$i; for(i=1; i<=NF; i++) if(tally[i]==0) tally[i]=1; for(i=1; i<=NF; i++) printf allele1[i]/tally[i]"\t"; printf "\n"; for(i=1; i<=NF; i++) printf allele2[i]/tally[i]"\t"; printf "\n"; for(i=1; i<=NF; i++) tally[i]=0}' MyData | sed 's/\t$//g'

Mas não sei como traduzi-lo para perl ou python >

    

Answer 1

Aqui está uma solução básica em Python (sem pacotes extravagantes de terceiros) que deve ser bem parecida com o que você está procurando:

#!/usr/bin/env python2

# snp.py

import sys

# Get the name of the data file from the command-line
data_file = sys.argv[1]

# Read and parse the data from the text file
data = []
with open(data_file, 'r') as file_handle:
    for line in file_handle:
        data.append([float(n) for n in line.split()])

# Get the number of rows and columns
rows = len(data)
cols = len(data[0])

# Iterate over adjacent pairs of rows
for r in range(rows//2):

    # Iterate over columns
    for c in range(cols):

        # Compute the sum of the two matching entries the pair of rows
        t = data[2*r][c] + data[2*r+1][c]
        if t:

            # Divide each entry by the sum of the pair
            data[2*r][c] /= t
            data[2*r+1][c] /= t

# Convert the data array back into formatted strings and print the results
for row in data:
    print(' '.join(['{0: <8}'.format(round(x,6)) for x in row]))

Isso funciona para você? Deve ser fácil ajustar se estiver um pouco fora da formatação.