Skip to content

Commit

Permalink
Add support for only calculating coverage for cut sites instead of th…
Browse files Browse the repository at this point in the history
…e whole fragment.

Add support for only calculating coverage for cut sites instead of the whole fragment.
  • Loading branch information
ghuls committed Jan 9, 2024
1 parent 044a19f commit b9a97d5
Show file tree
Hide file tree
Showing 3 changed files with 26 additions and 8 deletions.
3 changes: 2 additions & 1 deletion src/scatac_fragment_tools/cli/commands.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,6 +32,7 @@ def command_fragment_to_bigwigs(args):
fragments_df,
chrom_sizes,
args.bigwig_filename,
args.normalize == "yes",
args.normalize,
args.scaling_factor,
args.cut_sites,
)
22 changes: 16 additions & 6 deletions src/scatac_fragment_tools/cli/main.py
Original file line number Diff line number Diff line change
Expand Up @@ -55,12 +55,11 @@ def add_fragments_to_bigwig_subparser(
"-n",
"--normalize",
dest="normalize",
action="store",
type=str,
action="store_true",
required=False,
choices={"yes", "no"},
default="yes",
help='Normalize genome coverage data. Default: "yes".',
default=False,
help="Normalize genome coverage data by dividing by sequencing depth "
"(number of fragments) multiplied by 1 million.",
)
optional_arguments.add_argument(
"-s",
Expand All @@ -70,7 +69,18 @@ def add_fragments_to_bigwig_subparser(
type=float,
required=False,
default=1.0,
help='Scaling factor for genome coverage data. Default: "1.0".',
help="Scaling factor for genome coverage data. "
'If normalization is enabled, scaling is applied afterwards. Default: "1.0".',
)
optional_arguments.add_argument(
"-x",
"--cut-sites",
dest="cut_sites",
action="store_true",
required=False,
default=False,
help="Use 1 bp Tn5 cut sites (start and end of each fragment) instead of whole "
"fragment length for coverage calculation.",
)
optional_arguments.add_argument(
"--chrom-prefix",
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -7,13 +7,14 @@
import numpy as np
import polars as pl
import pyarrow as pa
import pyarrow.csv
import pyBigWig


def get_chromosome_sizes(chrom_sizes_filename: str):
chrom_sizes = {}

with open(chrom_sizes_filename, "r") as fh:
with open(chrom_sizes_filename, "r") as fh: # noqa: UP015
for line in fh:
line = line.rstrip("\n")

Expand Down Expand Up @@ -332,6 +333,12 @@ def fragments_to_bw(
starts, ends = (
per_chrom_fragments_dfs[chrom].select(["Start", "End"]).to_numpy().T
)
if cut_sites:
# Create cut site positions (for both start and end of a fragment).
starts, ends = (
np.hstack((starts, ends - 1)),
np.hstack((starts + 1, ends)),
)
chrom_arrays[chrom] = calculate_depth(chrom_sizes[chrom], starts, ends)
no_fragments += per_chrom_fragments_dfs[chrom].height

Expand Down

0 comments on commit b9a97d5

Please sign in to comment.