awkの沼につかった件について

2021年9月19日 04:05

前書き

私はプログラマーでもないし、ましてやPCなんて全く詳しくないのだが、時と場合によってはWindowsでしないような操作を要求されることがある。現在私が行っている研究テーマでは、Linuxとやらを使ったCUIで解析やらをしなくてはならない。Linuxといっても、GCPのVMインスタンスを使っているのでOSはDebianだ。Debianがどういうものかも知らない人間が適当にコマンドを打ち込んで解析なんてしてるのだから、我ながら恥ずかしい限りである。そんな人間が今回、awkの沼につかったことについて書いていく。この記事が読者の皆様の参考になれば幸いだ。（到底なるとは思えないが）

awkコマンドについて

awkはテキストファイルを行ごとに処理できるコマンドらしい。列をしていすることができ、何列目だけ表示したいとか、とある列ともう一つの列の数字を足したりなど様々なことができる。

1行目だけ表示させたいとき

awk '{print $1}' sample.txt

1行目と2行目を足して出力させたいとき

awk '{print $1 + $2}' sample.txt

さて、こんなコマンドを使って何がしたいのかと思われるかもしれない。実は私、生物系を専門とする大学に通っていて、そこで生命現象について解析を行っている。そこで使うパラメーターの一つが”TPM”であった。

TPMについて

”TPM は、サンプル中に全転写産物が 100 万個存在するときに、各転写産物に何個あたりの転写産物が存在するのかを表す値である。”

詳しい説明は省くが、要はこれでどれだけその遺伝子が発現してるか分かるらしい。なのでこのTPMを計算するために、csvファイルに記載されているリードカウントやら転写産物の長さやらを使ってawkを使おうと思ったのだ。

まずTt（転写産物tの1000bp当たりのリード数）を計算する。これは簡単で、リードカウントを長さで割って1000をかければよい。ファイルの構造は1列目に転写産物名、2列目長さ、3,4,5列目にリードカウントがあるため、まず3,4,5列目を合計して6列目に追加する。ここでpasteを使う。

まず合計を6列目に追加

paste TPM3.txt <(awk '{print $3 + $4 + $5}' TPM3.txt) > TPM4.txt

そして6列目を使って7列目にTtを追加

paste TPM4.txt <(awk '{print $6 * 1000 / $2}' TPM4.txt) > TPM5.txt

このファイル名（TPM3、TPM4）は自分で勝手につけた名前なのでお気になさらず。この番号のつけかた何号機って感じがして好きです。そして最後にTPMを計算しようとしたのだが、ここで苦労する。TPMの計算には全Ttの合計値を使うのだが、この計算にもawkが使える。

1列目の合計を出すとき

awk '{sum+=$1} END {print sum}' sample.txt

ならawkの中にawkを使えたら勝ちなのでは？と思って最初はこのようなことをしようとした。

paste TPM5.txt <(awk '{print $7 * 1000000 / <(awk '{sum+=$7} END {print sum}')}' TPM5.txt) > TPM6.txt

が、ダメ。どうにもawkの中にawkは使えなさそうということでSUMの値を変数にぶち込んで入れてみてはどうかとやってみた。

SUM=`awk ‘{sum+=$7} END {print sum}’ TPM5.txt` | paste TPM5.txt <(awk ‘{print $7 * 1000000 / $SUM }’ TPM5.txt ) > TPM6.txt

これもダメ。なんかそもそもSUMの値は数値として認識されてなさそう？ということでdeclare -i をしたり、計算結果から累乗で表示されてたためそうでない小数点表示にしたりなど色々やったが全てエラー。もう仕方なくSUMを小数点表示にしてから9行目に出力して、それを使って10行目にTPMを出力するという全くスマートじゃない方法を使おうかとまで思った。

phpの登場

そこで助っ人が登場、「phpとかpython使えばいいじゃない」と。

？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？

名前は聞いたことあるけど…何も分からないものをどう使えと？と言ったらコードを書いてくれました。

<?php
	// ファイル名：後で引数に変更
	$file_name = "TPM5.txt";
	// 全部読み出して、行単位で配列に入れます
	$ret_array = file( $file_name );
//	print_r($ret_array);

	$count = count($ret_array);
	$summary = 0;

	// １行単位でループ
	for ($i = 0; $i < $count; $i ++)
	{
		// １行取り出します
		$current = $ret_array[$i];
		// タブで分解
		$splits = preg_split("/\t/", $current);
		// 浮動小数点変換
		$t = (float)$splits[6];
		// 足し算
		$summary += $t;
	}

	// 保存用の配列を準備
	$result = [];

	// ２回目
	for ($i = 0; $i < $count; $i ++)
	{
		// １行取り出します
		$current = $ret_array[$i];
		// タブで分解
		$splits = preg_split("/\t/", $current);
		// 浮動小数点変換
		$t = (float)$splits[6];
		$splits[6] = (float)$splits[6];
		// tmpは１個分の計算値
		$tmp = $t / $summary * 1000000;
		// 配列に追加：７個目
		array_push($splits, $tmp);

		// 出力配列にまるごと追加
		array_push($result, $splits);
	}

//	print_r($result);
/*
	// 
	$fp = fopen("TPM6.txt", "w");
	for ($i = 0; $i < $result; $i ++)
	{
		$current = $ret_array[$i];
		for ($j = 0; $j < count($current); $j ++)
		{
			fwrite($fp, $current);
		}
	}
*/

	$fp = fopen('TPM6.txt', 'w');
	for ($i = 0; $i < $result; $i ++)
	{
		$current = $result[$i];
		fputcsv($fp, $current, "\t");
	}
	fclose($fp);

?>

ご丁寧に説明まで書いてくれた。やってることは同じだが、こういう計算はphpのほうが得意なのだろう。私も使い方を覚えなくてはならない…。

awkの問題点

・awkのなかにawkは使えない

・変数が使えない

これらは私の圧倒的技術不足の可能性も多いため、もし有識者がいらっしゃいましたらコメントしていただけると大変ありがたいです。

参考にした記事

この記事が気に入ったらサポートをしてみませんか？