我需要解析使用不同参数运行的化学程序的输出,并以特定格式组合感兴趣的信息。
该程序的每个输出文件如下表所示,它给出了特定残留物(此处为pH=0)上质子化和非质子化物种(残留物)的数量:
Residue Number State 0 State 1 State 2 State 3 State 4
-----------------------------------------------------------------------------------
Residue: GL4 7 0.000410 (0) 0.453512 (1) 0.004275 (1) 0.535908 (1) 0.005895 (1)
Residue: HIP 15 0.900000 (2) 0.080000 (1) 0.020000 (1)
Residue: AS4 18 0.010085 (0) 0.486042 (1) 0.004335 (1) 0.495922 (1) 0.003615 (1)
Residue: GL4 35 0.000000 (0) 0.581343 (1) 0.000360 (1) 0.368002 (1) 0.050295 (1)
Residue: AS4 48 0.022640 (0) 0.520073 (1) 0.018440 (1) 0.425152 (1) 0.013695 (1)
Residue: AS4 52 0.038725 (0) 0.517533 (1) 0.113676 (1) 0.280601 (1) 0.049465 (1)
Residue: AS4 66 1.000000 (0) 0.000000 (1) 0.000000 (1) 0.000000 (1) 0.000000 (1)
Residue: AS4 87 0.004295 (0) 0.439747 (1) 0.010535 (1) 0.524678 (1) 0.020745 (1)
Residue: AS4 101 0.000105 (0) 0.504673 (1) 0.013110 (1) 0.478517 (1) 0.003595 (1)
Residue: AS4 119 0.014240 (0) 0.488767 (1) 0.007100 (1) 0.483272 (1) 0.006620 (1)对于每个pH,我都有一个类似这样的文件(所有文件都有完全相同的残留物和状态,只是种群发生了变化)。现在我想提取所有残基的去质子化部分。去质子化分数对应于其编号后有(0)的群体:例如,在pH=0的GL4 7的情况下,它是0.000410 (对应于状态0),对于AS4 66,它是1.00000。实际上,除HIP 15外,所有残基的状态都是0.080000 :在这种情况下,去质子化分数用(1)表示,对应于状态1和2。在上面的例子中,它是0.080000+ 0.020000 = 0.1.
然后,我需要将这些来自不同文件的信息组合成一个文件,如下所示:
# pH GLU7 HIS15 ASP18 GLU35 ASP48 ASP52 ASP66 ASP87 ASP101 ASP119
0.000 0.000 0.100 0.010 0.000 0.023 0.039 1.000 0.004 0.000 0.014
1.000 0.006 0.140 0.098 0.000 0.276 0.312 1.000 0.015 0.002 0.069每列对应一个残基,每行对应一个pH (即来自单个文件的信息,这里我只显示来自两个文件的信息)。
我试着想出一些awk的一行,但我是一个初学者,我不确定如何继续。实际上,我不知道awk是不是最适合这项工作的工具。也许sed和grep或python会更好。我将需要使用许多不同的输出(但所有输出看起来都是一样的,尽管残差会发生变化)多次进行这种解析,所以我希望有一种方法可以使这一过程自动化,但具有一定的灵活性。
如果您有任何建议或意见,请不要犹豫,如果您能帮助我解决这个问题,我将不胜感激。
首先要感谢大家!
发布于 2016-12-04 07:49:10
您可以使用for循环将所有文件分类到一个文件中,并使用Stackoverflow中的前一个解决方案将行转置为列。
发布于 2016-12-04 07:53:46
不完全清楚您想要什么,但python的拆分函数可能会对您有用。如果在没有任何参数的情况下调用,它将根据空格进行拆分(将多个空格整理成一个)
以这行为例,
Residue: GL4 7 0.000410 (0) 0.453512 (1) 0.004275 (1) 0.535908 (1) 0.005895 (1)可以像这样分割,
a = 'Residue: GL4 7 0.000410 (0) 0.453512 (1) 0.004275 (1) 0.535908 (1) 0.005895 (1)'
l = a.split()
print l
['Residue:', 'GL4', '7', '0.000410', '(0)', '0.453512', '(1)', '0.004275', '(1)', '0.535908', '(1)', '0.005895', '(1)']然后,您可以访问所需的值并对其进行处理。对字符串调用float和int (例如,float('0.00410')应该会为你把它们转换成数字。对于'(1)',你可以做int('(1)'1:-1)
发布于 2016-12-05 18:42:34
这个awk脚本应该可以让您开始使用它。为了获得所需的输出,您必须将文件名替换为相应的pH值。我省略了不包含零状态的行,因为您没有指定如何处理这些行。
/^ Residue/ || /^-----/ { next; }
{
filenames[FILENAME] = 1;
columns[$2 " " $3] = 1;
for (i = 5; i <= NF; i = i + 2) {
if ($i == "(0)") {
data[$2 " " $3, FILENAME] = $(i-1);
}
}
}
END {
printf("%10s", "filename");
for (col in columns) {
printf("%10s", col);
}
print "";
for (filename in filenames) {
printf("%10s", filename);
for (col in columns) {
printf("%10s", data[col, filename]);
}
print "";
}
}https://stackoverflow.com/questions/40953749
复制相似问题