-
Notifications
You must be signed in to change notification settings - Fork 1
Expand file tree
/
Copy pathsync_cut.py
More file actions
104 lines (78 loc) · 3.06 KB
/
Copy pathsync_cut.py
File metadata and controls
104 lines (78 loc) · 3.06 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
from pathlib import Path
import argparse
def process_sync_file(input_sync, cut_sync, output_prefix):
"""
1. Cuts the input sync file to the first 70 SNP columns.
2. Splits the resulting file into one file per chromosome.
Parameters
----------
input_sync : str or Path
Path to the original input .sync file.
cut_sync : str or Path
Path where the cut .sync file will be written.
output_prefix : str or Path
Prefix for the chromosome-specific output files.
The chromosome name and '.sync' are appended automatically.
"""
input_sync = Path(input_sync)
cut_sync = Path(cut_sync)
output_prefix = Path(output_prefix)
# ---------------------------------------------------------
# Step 1: Cut sync file to the first 70 SNP columns
# ---------------------------------------------------------
print(f"Cutting sync file:\n {input_sync}\n-> {cut_sync}")
with input_sync.open("r") as infile, cut_sync.open("w") as outfile:
for line in infile:
columns = line.rstrip("\n").split("\t")
# Keep chromosome, position, reference base
# and the first 70 SNP columns
outfile.write("\t".join(columns[:73]) + "\n")
print("Cutting completed.\n")
# ---------------------------------------------------------
# Step 2: Split cut sync file by chromosome
# ---------------------------------------------------------
print(f"Splitting sync file by chromosome:\n {cut_sync}")
current_chromosome = None
output_file = None
try:
with cut_sync.open("r") as infile:
for line in infile:
columns = line.rstrip("\n").split("\t")
chromosome = columns[0]
# Open a new output file when chromosome changes
if chromosome != current_chromosome:
if output_file is not None:
output_file.close()
current_chromosome = chromosome
output_path = Path(
f"{output_prefix}{current_chromosome}.sync"
)
print(f" Writing chromosome {current_chromosome} -> {output_path}")
output_file = output_path.open("w")
output_file.write(line)
finally:
if output_file is not None:
output_file.close()
print("\nDone.")
if __name__ == "__main__":
parser = argparse.ArgumentParser(
description="Cut a sync file to 70 columns and split it by chromosome."
)
parser.add_argument(
"input_sync",
help="Path to the original input .sync file."
)
parser.add_argument(
"cut_sync",
help="Path to the output .sync file containing the first 70 SNP columns."
)
parser.add_argument(
"output_prefix",
help="Prefix for chromosome-specific output files."
)
args = parser.parse_args()
process_sync_file(
input_sync=args.input_sync,
cut_sync=args.cut_sync,
output_prefix=args.output_prefix,
)