Repository: akikuno/DAJIN2 (language: python)

Issue:
Test for Dorado basecaller
## Describe the features you want.

Current DAJIN2 assumes base calling with Guppy, but it is necessary to confirm whether it works similarly with Dorado.

## Describe the solution you'd like if you have any.

Try basecall and demultiplex using [Dorado](https://github.com/nanoporetech/dorado?tab=readme-ov-file)

## DAJIN2 version

0.4.6

## Additional context

None


Discussion / hints:
In Dorado, unlike Guppy, the output is in BAM format.
![image](https://github.com/akikuno/DAJIN2/assets/15861316/cc94ee7e-a5d9-4b88-8aba-d9dbdcbe088a)

Considering that Dorado is likely to become mainstream in the future, it is better to primarily support Dorado. **The input should accept BAM files and convert them to FASTQ internally.**

For Guppy, the FASTQ files are output inside the `barcodeXX` directory.
**In the case of Guppy, it might be good to allow directory input by adding a `--guppy`option.**

![image](https://github.com/akikuno/DAJIN2/assets/15861316/b842138b-e6b7-4fbb-a708-019a819961d4)


Fix patch:
diff --git a/README.md b/README.md
index ddff025..7f12195 100644
--- a/README.md
+++ b/README.md
@@ -66,13 +66,12 @@ pip install DAJIN2
 
 ### Required Files
 
-#### FASTQ Files for Sample and Control
+#### FASTQ/FASTA/BAM Files for Sample and Control
 
-In DAJIN2, a **control that has not undergone genome editing** is necessary to detect genome-editing-specific mutations. Specify a directory containing the FASTQ files (both gzip compressed and uncompressed) of the genome editing sample and control.
-
-<!-- [Nanopore Guppy](https://community.nanoporetech.com/docs/prepare/library_prep_protocols/Guppy-protocol) -->
-After base calling with Guppy, the following file structure will be output:
+In DAJIN2, a **control that has not undergone genome editing** is necessary to detect genome-editing-specific mutations. Specify a directory containing the FASTQ/FASTA (both gzip compressed and uncompressed) or BAM files of the genome editing sample and control.
 
+##### [Basecalling with Guppy](https://community.nanoporetech.com/docs/prepare/library_prep_protocols/Guppy-protocol/v/gpb_2003_v1_revax_14dec2018/guppy-software-overview)
+After basecalling with Guppy, the following file structure will be output:
 
 ```text
 fastq_pass
@@ -86,12 +85,49 @@ fastq_pass
     └── fastq_runid_b347657c88dced2d15bf90ee6a1112a3ae91c1af_11_0.fastq.gz
 ```
 
-Assuming barcode01 as the control and barcode02 as the sample, specify each directory as follows:
+Assuming barcode01 is the control and barcode02 is the sample, the respective directories are specified as follows:
+
++ Control: `fastq_pass/barcode01`
++ Sample: `fastq_pass/barcode02`
+
+
+##### [Basecalling with Dorado](https://github.com/nanoporetech/dorado)
+
+For basecalling with Dorado ([`dorado demux`](https://github.com/nanoporetech/dorado?tab=readme-ov-file#barcode-classification)), the following file structure will be output:
+
+```text
+dorado_demultiplex
+├── EXP-PBC096_barcode01.bam
+└── EXP-PBC096_barcode02.bam
+```
+
+> [!IMPORTANT]
+> Store each BAM file in a separate directory. The directory names can be set arbitrarily.
+
+```text
+dorado_demultiplex
+├── barcode01
+│   └── EXP-PBC096_barcode01.bam
+└── barcode02
+    └── EXP-PBC096_barcode02.bam
+```
+
+Similarly, store the FASTA files outputted after sequence error correction with [`dorado correct`](https://github.com/nanoporetech/dorado) in separate directories.
+
+```text
+dorado_correct
+├── barcode01
+│   └── EXP-PBC096_barcode01.fasta
+└── barcode02
+    └── EXP-PBC096_barcode02.fasta
+```
+
+Assuming barcode01 is the control and barcode02 is the sample, the respective directories are specified as follows:
 
-+ Control: fastq_pass/barcode01
-+ Sample: fastq_pass/barcode01
++ Control: `dorado_demultiplex/barcode01` / `dorado_correct/barcode01`
++ Sample: `dorado_demultiplex/barcode02` / `dorado_correct/barcode02`
 
-#### FASTA File for Allele
+#### FASTA File Including Anticipated Allele Sequences
 
 The FASTA file should contain descriptions of the alleles anticipated as a result of genome editing.
 
@@ -123,8 +159,8 @@ DAJIN2 <-s|--sample> <-c|--control> <-a|--allele> <-n|--name> \
   [-g|--genome] [-t|--threads] [-h|--help] [-v|--version]
 
 Options:
--s, --sample              Specify the path to the directory containing sample FASTQ files.
--c, --control             Specify the path to the directory containing control FASTQ files.
+-s, --sample              Specify the path to the directory containing sample FASTQ/FASTA/BAM files.
+-c, --control             Specify the path to the directory containing control FASTQ/FASTA/BAM files.
 -a, --allele              Specify the path to the FASTA file.
 -n, --name (Optional)     Set the output directory name. Default: 'Results'.
 -g, --genome (Optional)   Specify the reference UCSC genome ID (e.g., hg38, mm39). Default: '' (empty string).
@@ -181,7 +217,7 @@ DAJIN2 \
 
 ### Batch Processing
 
-By using the `batch` subcommand, you can process multiple FASTQ files simultaneously.  
+By using the `batch` subcommand, you can process multiple files simultaneously.  
 For this purpose, a CSV or Excel file consolidating the sample information is required.  
 
 > [!NOTE]
diff --git a/docs/README_JP.md b/docs/README_JP.md
index 68083eb..32c461a 100644
--- a/docs/README_JP.md
+++ b/docs/README_JP.md
@@ -40,7 +40,7 @@ conda activate env-dajin2
 > 現状、[BiocondaがPython 3.11以上に対応していない](https://github.com/bioconda/bioconda-recipes/issues/37805)ため、DAJIN2はPython 3.8 から 3.10までをサポートしています。
 
 > [!NOTE]
-> [Apple SiliconはBiocondaチャンネルに対応していない](https://github.com/bioconda/bioconda-recipes/issues/37068#issuecomment-1257790919)ため、以下のようにRoseeta2経由でインストールをしてください
+> 現状、[BiocondaがApple Siliconに対応していない](https://github.com/bioconda/bioconda-recipes/issues/37068#issuecomment-1257790919)ため、以下のようにRoseeta2経由でインストールを行ってください
 > ```bash
 > CONDA_SUBDIR=osx-64 conda create -n env-dajin2 -c conda-forge -c bioconda python=3.10 DAJIN2 -y
 > conda activate env-dajin2
@@ -62,14 +62,14 @@ pip install DAJIN2
 
 ### 必要なファイル
 
-#### サンプルおよびコントロールのFASTQファイル
+#### サンプルおよびコントロールファイル
 
-DAJIN2では、ゲノム編集特異的な変異を検出するために、**ゲノム編集を受けていないコントロール**が必要です。  
-ゲノム編集サンプルとコントロールのFASTQファイル（Gzip圧縮・非圧縮どちらも対応可能）を含むディレクトリを指定します。
+DAJIN2では、ゲノム編集特異的な変異を検出するために、**ゲノム編集を受けていないコントロールサンプル**が必要です。  
+ゲノム編集サンプルとコントロールのFASTQ/FASTA（gzip圧縮・非圧縮どちらも対応可能）、またはBAMファイルを含むディレクトリを指定します。
 
-<!-- [Nanopore Guppy](https://community.nanoporetech.com/docs/prepare/library_prep_protocols/Guppy-protocol) -->
-Guppyによるベースコール処理後、以下のようなファイル構成が出力されます：
+##### [Guppy](https://community.nanoporetech.com/docs/prepare/library_prep_protocols/Guppy-protocol/v/gpb_2003_v1_revax_14dec2018/guppy-software-overview)によるベースコール
 
+Guppyによるベースコール後、以下のようなファイル構成が出力されます：
 
 ```text
 fastq_pass
@@ -86,9 +86,46 @@ fastq_pass
 上記のbarcode01をコントロール、barcode02をサンプルと仮定すると、それぞれのディレクトリは下記の通りに指定します：
 
 + コントロール: `fastq_pass/barcode01`
-+ サンプル: `fastq_pass/barcode01`
++ サンプル: `fastq_pass/barcode02`
 
-#### FASTAファイル
+##### [Dorado](https://github.com/nanoporetech/dorado)によるベースコール
+
+Doradoによるベースコール（[`dorado demux`](https://github.com/nanoporetech/dorado?tab=readme-ov-file#barcode-classification)）においては、BAMファイルが出力されます：
+
+```text
+dorado_demultiplex
+├── EXP-PBC096_barcode01.bam
+└── EXP-PBC096_barcode02.bam
+```
+
+> [!IMPORTANT]
+> 各bamファイルを別々のディレクトリに格納してください。ディレクトリ名は任意です。
+
+```text
+dorado_demultiplex
+├── barcode01
+│   └── EXP-PBC096_barcode01.bam
+└── barcode02
+    └── EXP-PBC096_barcode02.bam
+```
+
+[`dorado correct`](https://github.com/nanoporetech/dorado)によるシークエンスエラー補正後に出力されるFASTAファイルも同様に、別々のディレクトリに格納してください。
+
+```text
+dorado_correct
+├── barcode01
+│   └── EXP-PBC096_barcode01.fasta
+└── barcode02
+    └── EXP-PBC096_barcode02.fasta
+```
+
+上記のbarcode01をコントロール、barcode02をサンプルと仮定すると、それぞれのディレクトリは下記の通りに指定します：
+
++ コントロール: `dorado_demultiplex/barcode01` / `dorado_correct/barcode01`
++ サンプル: `dorado_demultiplex/barcode02` / `dorado_correct/barcode02`
+
+
+#### 想定されるアレル配列を含むFASTAファイル
 
 FASTAファイルには、ゲノム編集によって想定されるアレルを記述します。
 
@@ -120,8 +157,8 @@ DAJIN2 <-c|--control> <-s|--sample> <-a|--allele> <-n|--name> \
   [-g|--genome] [-t|--threads] [-h|--help] [-v|--version]
 
 引数:
-  -c, --control             コントロールのFASTQファイルが格納されたディレクトリのパス
-  -s, --sample              サンプルのFASTQファイルが格納されたディレクトリのパス
+  -c, --control             コントロールのFASTQ/FASTA/BAMファイルが格納されたディレクトリのパス
+  -s, --sample              サンプルのFASTQ/FASTA/BAMファイルが格納されたディレクトリのパス
   -a, --allele              ゲノム編集によって想定されるアレルを記載したFASTAファイルのパス
   -n, --name (オプション)     出力ディレクトリの名前 [デフォルト: Results]
   -g, --genome (オプション)   参照ゲノムのID (e.g hg38, mm39) [デフォルト: '']
diff --git a/docs/RELEASE.md b/docs/RELEASE.md
index 2518005..a404ff0 100644
--- a/docs/RELEASE.md
+++ b/docs/RELEASE.md
@@ -11,16 +11,46 @@
 
 <!-- 💡 ToDo
 
-- VCF、PDFを出力する
-- 逆位アレルでの検証を加える
-- nCATSがほしい…
-- Flaskではなく、streamlitでGUIを作る
-
  -->
 
 <!-- ############################################################# # -->
 
-# v0.5.0 (2024-MM-DD)
+## 💥 Breaking
+
++ Enable to accept additional file formats as an input. Issue: #37
+  + FASTA [[Commit Detail](https://github.com/akikuno/DAJIN2/commit/ee6d392cd51649c928bd604acafbab4b9d28feb1)]
+  + BAM [[Commit Detail](https://github.com/akikuno/DAJIN2/commit/xxx)]
+
+## 📝 Documentation
+
++ Add a description of the procedure for accepting files generated by Dorado basecaller as input. Issue: #37 [[Commit Detail](https://github.com/akikuno/DAJIN2/commit/c9ebc020fa60980ba7aaaf9295975775ec07da6d)]
+
+
+## 🔧 Maintenance
+
++ Specify the Python version to be between 3.8 and 3.10. [[Commit Detail](https://github.com/akikuno/DAJIN2/commit/5fae947eff7da0f7e1ed5e4ff3f95c911fd9f646)]
+
++ Change `mutation_exporter.report_mutations` to return list[list[str]]. Update the tests accordingly. [[Commit Detail](https://github.com/akikuno/DAJIN2/commit/7153cb143d621e136ca94bfe6b391f1d7b61d438)]
+
++ Apply formatting with Ruff [[Commit Detail](https://github.com/akikuno/DAJIN2/commit/xxaec9b697863ef06b4e86e248bebde6616f4eb54ex)]
+
+## 🐛 Bug Fixes
+
++ Add `reallocate_insertion_within_deletion` into `report.mutation_exporter` and reflected it in the mutation info. [[Commit Detail](https://github.com/akikuno/DAJIN2/commit/ed6a96e01bb40c77df9cd3a17a4c29524684b6f1)]
+
+
+<!-- ############################################################# # -->
+
+
+
+-------------------------------------------------------------
+
+# Past Releases
+
+<!--  ------------------------------------------------------------- -->
+
+<details>
+<summary> v0.5.0 (2024-06-05) </summary>
 
 ## 📝 Documentation
 
@@ -84,17 +114,7 @@ Update the algorithms of `cssplits_hander.reallocate_insertion_within_deletion`
 + Modified the system to separate intermediate files using a directory structure instead of underscores (`_`), ensuring that no errors occur even if users use allele names containing underscores [[Commit Detail](https://github.com/akikuno/DAJIN2/commit/f70948315a114e0c182895ba4320233f26fc1025)]
   + Thank you @geedrn for reporting the issue #39!
 
-
-
-<!-- ############################################################# # -->
-
-
-
--------------------------------------------------------------
-
-# Past Releases
-
-<!--  ------------------------------------------------------------- -->
+</details>
 
 <details>
 <summary> v0.4.6 (2024-05-17) </summary>
diff --git a/pyproject.toml b/pyproject.toml
index 1afb3ad..2654223 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -4,7 +4,7 @@ build-backend = "poetry.core.masonry.api"
 
 [tool.poetry]
 name = "DAJIN2"
-version = "0.5.0"
+version = "0.5.1"
 description = "One-step genotyping tools for targeted long-read sequencing"
 authors = ["Akihiro Kuno <akuno@md.tsukuba.ac.jp>"]
 readme = "README.md"
@@ -29,7 +29,7 @@ include = [
 ]
 
 [tool.poetry.dependencies]
-python = "^3.8"
+python = ">=3.8, <3.11"
 numpy = ">=1.24.0"
 scipy = ">=1.10.0"
 pandas = ">=1.0.0"
@@ -51,3 +51,41 @@ ruptures = ">=1.1.8"
 [tool.poetry.scripts]
 DAJIN2 = "DAJIN2.main:execute"
 
+
+[tool.ruff]
+line-length = 119
+
+[tool.ruff.lint]
+select = [
+    "E",  # pycodestyle errors
+    "W",  # pycodestyle warnings
+    "F",  # pyflakes
+    "I",  # isort
+    "B",  # flake8-bugbear
+    "C4",  # flake8-comprehensions
+    "UP",  # pyupgrade
+]
+ignore = [
+    "E501",  # line too long, handled by black
+    "B008",  # do not perform function calls in argument defaults
+    "C901",  # too complex
+    "W191",  # indentation contains tabs
+    "B904", # ignore errors for raise ... from ... not being used
+]
+
+[tool.ruff.lint.per-file-ignores]
+"__init__.py" = ["F401"]
+
+[tool.ruff.lint.isort]
+combine-as-imports = true
+known-first-party = ["musubi_restapi"]
+section-order = ["future", "standard-library", "third-party", "first-party", "local-folder"]
+split-on-trailing-comma = true
+
+[tool.ruff.format]
+quote-style = "double"
+
+
+[tool.ruff.lint.pyupgrade]
+# Settings for Python 3.8 compatibility
+keep-runtime-typing = true
diff --git a/src/DAJIN2/core/classification/allele_merger.py b/src/DAJIN2/core/classification/allele_merger.py
index 00b812d..835128e 100644
--- a/src/DAJIN2/core/classification/allele_merger.py
+++ b/src/DAJIN2/core/classification/allele_merger.py
@@ -1,8 +1,7 @@
 from __future__ import annotations
 
-from itertools import groupby
 from collections import defaultdict
-
+from itertools import groupby
 
 ##########################################################
 # merge minor alleles
diff --git a/src/DAJIN2/core/classification/classifier.py b/src/DAJIN2/core/classification/classifier.py
index a3074d1..dc618c1 100644
--- a/src/DAJIN2/core/classification/classifier.py
+++ b/src/DAJIN2/core/classification/classifier.py
@@ -1,10 +1,10 @@
 from __future__ import annotations
 
-from pathlib import Path
 from itertools import groupby
+from pathlib import Path
 
-from DAJIN2.utils import io
 from DAJIN2.core.classification.allele_merger import merge_minor_alleles
+from DAJIN2.utils import io
 
 
 def calc_match(cssplit: str) -> float:
diff --git a/src/DAJIN2/core/clustering/__init__.py b/src/DAJIN2/core/clustering/__init__.py
index afaabe8..4ff1789 100644
--- a/src/DAJIN2/core/clustering/__init__.py
+++ b/src/DAJIN2/core/clustering/__init__.py
@@ -1,3 +1,3 @@
+from DAJIN2.core.clustering.appender import add_labels, add_percent, add_readnum
 from DAJIN2.core.clustering.label_extractor import extract_labels
-from DAJIN2.core.clustering.appender import add_labels, add_readnum, add_percent
 from DAJIN2.core.clustering.label_updator import update_labels
diff --git a/src/DAJIN2/core/clustering/clustering.py b/src/DAJIN2/core/clustering/clustering.py
index c50646c..a336336 100644
--- a/src/DAJIN2/core/clustering/clustering.py
+++ b/src/DAJIN2/core/clustering/clustering.py
@@ -1,18 +1,18 @@
 from __future__ import annotations
 
-from pathlib import Path
-from itertools import chain
 from collections import Counter
+from itertools import chain
+from pathlib import Path
 
 import numpy as np
+from scipy.sparse import csr_matrix, spmatrix
 from sklearn import metrics
 from sklearn.cluster import BisectingKMeans
-from scipy.sparse import csr_matrix, spmatrix
 
-from DAJIN2.utils import io, config
 from DAJIN2.core.clustering.label_merger import merge_labels
 from DAJIN2.core.clustering.score_handler import subset_scores
 from DAJIN2.core.clustering.strand_bias_handler import remove_biased_clusters
+from DAJIN2.utils import config, io
 
 config.set_warnings_ignore()
 
diff --git a/src/DAJIN2/core/clustering/label_extractor.py b/src/DAJIN2/core/clustering/label_extractor.py
index d74482f..eb500f4 100644
--- a/src/DAJIN2/core/clustering/label_extractor.py
+++ b/src/DAJIN2/core/clustering/label_extractor.py
@@ -1,15 +1,14 @@
 from __future__ import annotations
 
 import uuid
-
-from pathlib import Path
 from itertools import groupby
+from pathlib import Path
 
-from DAJIN2.utils import io
-from DAJIN2.core.clustering.score_handler import make_score, annotate_score
+from DAJIN2.core.clustering.clustering import return_labels
 from DAJIN2.core.clustering.label_updator import relabel_with_consective_order
+from DAJIN2.core.clustering.score_handler import annotate_score, make_score
 from DAJIN2.core.clustering.strand_bias_handler import is_strand_bias
-from DAJIN2.core.clustering.clustering import return_labels
+from DAJIN2.utils import io
 
 
 def extract_labels(classif_sample, TEMPDIR, SAMPLE_NAME, CONTROL_NAME) -> list[dict[str]]:
diff --git a/src/DAJIN2/core/clustering/label_merger.py b/src/DAJIN2/core/clustering/label_merger.py
index 42cd17c..26bf8e4 100644
--- a/src/DAJIN2/core/clustering/label_merger.py
+++ b/src/DAJIN2/core/clustering/label_merger.py
@@ -1,6 +1,7 @@
 from __future__ import annotations
 
 from collections import Counter
+
 import numpy as np
 
 
diff --git a/src/DAJIN2/core/clustering/score_handler.py b/src/DAJIN2/core/clustering/score_handler.py
index db71036..6e0db9e 100644
--- a/src/DAJIN2/core/clustering/score_handler.py
+++ b/src/DAJIN2/core/clustering/score_handler.py
@@ -1,8 +1,8 @@
 from __future__ import annotations
 
-from typing import Generator
-from itertools import groupby
 from collections import Counter
+from itertools import groupby
+from typing import Generator
 
 from DAJIN2.core.clustering.kmer_generator import generate_mutation_kmers
 
diff --git a/src/DAJIN2/core/clustering/strand_bias_handler.py b/src/DAJIN2/core/clustering/strand_bias_handler.py
index f59832b..789d7fe 100644
--- a/src/DAJIN2/core/clustering/strand_bias_handler.py
+++ b/src/DAJIN2/core/clustering/strand_bias_handler.py
@@ -1,5 +1,13 @@
 from __future__ import annotations
 
+from collections import defaultdict
+from pathlib import Path
+from typing import Generator
+
+from sklearn.tree import DecisionTreeClassifier
+
+from DAJIN2.utils import io
+
 """
 Nanopore sequencing results often results in strand specific mutations even though the mutation is not strand specific, thus they are considered as sequencing errors and should be removed.
 
@@ -8,13 +16,6 @@
 Re-allocates reads belonging to clusters with strand bias to clusters without strand bias.
 """
 
-from pathlib import Path
-from typing import Generator
-from collections import defaultdict
-from sklearn.tree import DecisionTreeClassifier
-
-from DAJIN2.utils import io
-
 # Constants
 STRAND_BIAS_LOWER_LIMIT = 0.1
 STRAND_BIAS_UPPER_LIMIT = 0.9
diff --git a/src/DAJIN2/core/consensus/__init__.py b/src/DAJIN2/core/consensus/__init__.py
index 5a8496a..3f2d5b2 100644
--- a/src/DAJIN2/core/consensus/__init__.py
+++ b/src/DAJIN2/core/consensus/__init__.py
@@ -1,6 +1,4 @@
+from DAJIN2.core.consensus.clust_formatter import downsample_by_label, remove_minor_alleles
 from DAJIN2.core.consensus.consensus import call_consensus
-from DAJIN2.core.consensus.name_handler import call_allele_name
-from DAJIN2.core.consensus.name_handler import update_key_by_allele_name
-from DAJIN2.core.consensus.name_handler import add_key_by_allele_name
-from DAJIN2.core.consensus.clust_formatter import remove_minor_alleles, downsample_by_label
 from DAJIN2.core.consensus.mutation_extractor import cache_mutation_loci
+from DAJIN2.core.consensus.name_handler import add_key_by_allele_name, call_allele_name, update_key_by_allele_name
diff --git a/src/DAJIN2/core/consensus/clust_formatter.py b/src/DAJIN2/core/consensus/clust_formatter.py
index 5c7c619..47187cd 100644
--- a/src/DAJIN2/core/consensus/clust_formatter.py
+++ b/src/DAJIN2/core/consensus/clust_formatter.py
@@ -1,8 +1,8 @@
 from __future__ import annotations
 
 import random
-from itertools import groupby
 from collections import defaultdict
+from itertools import groupby
 
 
 def remove_minor_alleles(clust_sample: list[dict]) -> list[dict]:
diff --git a/src/DAJIN2/core/consensus/consensus.py b/src/DAJIN2/core/consensus/consensus.py
index 28eaba0..7ff65e7 100644
--- a/src/DAJIN2/core/consensus/consensus.py
+++ b/src/DAJIN2/core/consensus/consensus.py
@@ -1,14 +1,13 @@
 from __future__ import annotations
 
-from pathlib import Path
+from collections import defaultdict
 from dataclasses import dataclass
 from itertools import groupby
-from collections import defaultdict
+from pathlib import Path
 
 from DAJIN2.utils import io
 from DAJIN2.utils.cssplits_handler import call_sequence
 
-
 ###########################################################
 # call position weight matrix (cons_pergentage)
 ###########################################################
@@ -98,11 +97,10 @@ class ConsensusKey:
 
 
 def call_consensus(tempdir: Path, sample_name: str, clust_sample: list[dict]) -> tuple[dict[list], dict[str]]:
-
     clust_sample.sort(key=lambda x: [x["ALLELE"], x["LABEL"]])
 
-    cons_percentages = dict()
-    cons_sequences = dict()
+    cons_percentages = {}
+    cons_sequences = {}
 
     for (allele, label), group in groupby(clust_sample, key=lambda x: [x["ALLELE"], x["LABEL"]]):
         clust = list(group)
diff --git a/src/DAJIN2/core/consensus/mutation_extractor.py b/src/DAJIN2/core/consensus/mutation_extractor.py
index c81d9ad..033b767 100644
--- a/src/DAJIN2/core/consensus/mutation_extractor.py
+++ b/src/DAJIN2/core/consensus/mutation_extractor.py
@@ -1,14 +1,14 @@
 from __future__ import annotations
 
-from pathlib import Path
 from itertools import groupby
+from pathlib import Path
 
 import numpy as np
 from sklearn.cluster import MiniBatchKMeans
 
-from DAJIN2.utils import io
-from DAJIN2.core.preprocess.mutation_extractor import summarize_indels, extract_mutation_loci, minimize_mutation_counts
 from DAJIN2.core.consensus.similarity_searcher import cache_selected_control_by_similarity
+from DAJIN2.core.preprocess.mutation_extractor import extract_mutation_loci, minimize_mutation_counts, summarize_indels
+from DAJIN2.utils import io
 
 """
 Most of the code reuses `preprocess.cache_mutation_loci`.
@@ -25,7 +25,7 @@ def get_thresholds(path_indels_normalized_sample, path_indels_normalized_control
     indels_normalized_sample = io.load_pickle(path_indels_normalized_sample)
     indels_normalized_control = io.load_pickle(path_indels_normalized_control)
     indels_normalized_minimize_control = minimize_mutation_counts(indels_normalized_control, indels_normalized_sample)
-    thresholds = dict()
+    thresholds = {}
     for mut in {"+", "-", "*"}:
         values_sample = indels_normalized_sample[mut]
         values_control = indels_normalized_minimize_control[mut]
diff --git a/src/DAJIN2/core/consensus/name_handler.py b/src/DAJIN2/core/consensus/name_handler.py
index 2a54245..9dfec87 100644
--- a/src/DAJIN2/core/consensus/name_handler.py
+++ b/src/DAJIN2/core/consensus/name_handler.py
@@ -1,6 +1,7 @@
 from __future__ import annotations
 
 import re
+
 from DAJIN2.core.consensus.consensus import ConsensusKey
 
 
@@ -69,7 +70,7 @@ def call_allele_name(
 
 
 def update_key_by_allele_name(cons: dict, allele_names: dict[int, str]) -> dict:
-    cons_update = dict()
+    cons_update = {}
     for key in cons:
         old_allele = cons[key]
         new_allele = allele_names[key.label]
diff --git a/src/DAJIN2/core/consensus/similarity_searcher.py b/src/DAJIN2/core/consensus/similarity_searcher.py
index 698cbf6..5f3bee3 100644
--- a/src/DAJIN2/core/consensus/similarity_searcher.py
+++ b/src/DAJIN2/core/consensus/similarity_searcher.py
@@ -1,10 +1,9 @@
 from __future__ import annotations
 
-from pathlib import Path
 from collections import defaultdict
+from pathlib import Path
 
 import numpy as np
-
 from sklearn.neighbors import LocalOutlierFactor
 
 from DAJIN2.utils import io
@@ -28,7 +27,7 @@ def onehot_by_mutations(midsv_sample: list[dict]) -> dict[str, np.ndarray]:
 def calculate_percentage(
     mut_onehot_sample: dict[str, np.ndarray], coverage_match: np.ndarray[int]
 ) -> dict[str, np.ndarray]:
-    mut_percentage = dict()
+    mut_percentage = {}
     for mut, onehot in mut_onehot_sample.items():
         x = np.sum(onehot, axis=0) / coverage_match
         mut_percentage[mut] = np.where(np.isnan(x), 0, x)
@@ -36,14 +35,14 @@ def calculate_percentage(
 
 
 def get_values_to_mask(mut_percentage_sample: dict[str, np.ndarray], threshold=0.5) -> dict[str, np.ndarray[float]]:
-    mask = dict()
+    mask = {}
     for mut, percentage in mut_percentage_sample.items():
         mask[mut] = np.where(percentage > threshold, 0, percentage)
     return mask
 
 
 def apply_mask(mut_onehot: dict[str, np.ndarray], mask_sample: dict[str, np.ndarray[float]]):
-    mut_onehot_masked = dict()
+    mut_onehot_masked = {}
     for mut, onehot in mut_onehot.items():
         mut_onehot_masked[mut] = onehot * mask_sample[mut]
     return mut_onehot_masked
@@ -52,7 +51,7 @@ def apply_mask(mut_onehot: dict[str, np.ndarray], mask_sample: dict[str, np.ndar
 def identify_normal_reads(
     mut_onehot_sample_masked: dict[str, np.ndarray], mut_onehot_control_masked: dict[str, np.ndarray]
 ) -> list[bool]:
-    mutation_comparisons = dict()
+    mutation_comparisons = {}
     for mut in {"+", "-", "*"}:
         values_sample = mut_onehot_sample_masked[mut]
         values_control = mut_onehot_control_masked[mut]
diff --git a/src/DAJIN2/core/core.py b/src/DAJIN2/core/core.py
index 5f6dead..892f2f9 100644
--- a/src/DAJIN2/core/core.py
+++ b/src/DAJIN2/core/core.py
@@ -1,13 +1,12 @@
 from __future__ import annotations
 
-import shutil
 import logging
-
+import shutil
 from pathlib import Path
 
-from DAJIN2.utils import io, fastx_handler
 from DAJIN2.core import classification, clustering, consensus, preprocess, report
 from DAJIN2.core.preprocess.input_formatter import FormattedInputs
+from DAJIN2.utils import fastx_handler, io
 
 logger = logging.getLogger(__name__)
 
@@ -18,7 +17,6 @@
 
 
 def execute_control(arguments: dict):
-
     logger.info(f"{arguments['control']} is now processing...")
 
     ###########################################################
@@ -42,7 +40,7 @@ def execute_control(arguments: dict):
     # Merge fastq files
     ###########################################################
 
-    fastx_handler.save_concatenated_fastx(ARGS.tempdir, ARGS.path_control, ARGS.control_name)
+    fastx_handler.save_inputs_as_single_fastx(ARGS.tempdir, ARGS.path_control, ARGS.control_name)
     # Save subsetted fastq if the read number is too large (> 10,000 reads)
     fastx_handler.save_subsetted_fastx(
         Path(ARGS.tempdir, ARGS.control_name, "fastq", f"{ARGS.control_name}.fastq.gz"), num_reads=10_000
@@ -80,7 +78,7 @@ def execute_control(arguments: dict):
     ###########################################################
     # Finish call
     ###########################################################
-    logger.info(f"\N{teacup without handle} {arguments['control']} is finished!")
+    logger.info(f"\N{TEACUP WITHOUT HANDLE} {arguments['control']} is finished!")
 
 
 def execute_sample(arguments: dict):
@@ -100,7 +98,7 @@ def execute_sample(arguments: dict):
     # Merge fastq files
     # ============================================================
 
-    fastx_handler.save_concatenated_fastx(ARGS.tempdir, ARGS.path_sample, ARGS.sample_name)
+    fastx_handler.save_inputs_as_single_fastx(ARGS.tempdir, ARGS.path_sample, ARGS.sample_name)
 
     # ============================================================
     # Mapping with mappy
@@ -231,4 +229,4 @@ def execute_sample(arguments: dict):
     # Finish call
     ###########################################################
 
-    logger.info(f"\N{teacup without handle} {arguments['sample']} is finished!")
+    logger.info(f"\N{TEACUP WITHOUT HANDLE} {arguments['sample']} is finished!")
diff --git a/src/DAJIN2/core/preprocess/__init__.py b/src/DAJIN2/core/preprocess/__init__.py
index 73d7989..51d791f 100644
--- a/src/DAJIN2/core/preprocess/__init__.py
+++ b/src/DAJIN2/core/preprocess/__init__.py
@@ -1,9 +1,9 @@
-from DAJIN2.core.preprocess.cache_checker import exists_cached_hash, exists_cached_genome
-from DAJIN2.core.preprocess.genome_fetcher import fetch_coordinates, fetch_chromosome_size
-from DAJIN2.core.preprocess.mapping import generate_sam
-from DAJIN2.core.preprocess.directory_manager import create_temporal_directories, create_report_directories
+from DAJIN2.core.preprocess.cache_checker import exists_cached_genome, exists_cached_hash
+from DAJIN2.core.preprocess.directory_manager import create_report_directories, create_temporal_directories
+from DAJIN2.core.preprocess.genome_fetcher import fetch_chromosome_size, fetch_coordinates
 from DAJIN2.core.preprocess.input_formatter import format_inputs
-from DAJIN2.core.preprocess.midsv_caller import generate_midsv
+from DAJIN2.core.preprocess.insertions_to_fasta import generate_insertion_fasta
 from DAJIN2.core.preprocess.knockin_handler import extract_knockin_loci
+from DAJIN2.core.preprocess.mapping import generate_sam
+from DAJIN2.core.preprocess.midsv_caller import generate_midsv
 from DAJIN2.core.preprocess.mutation_extractor import cache_mutation_loci
-from DAJIN2.core.preprocess.insertions_to_fasta import generate_insertion_fasta
diff --git a/src/DAJIN2/core/preprocess/cache_checker.py b/src/DAJIN2/core/preprocess/cache_checker.py
index c6d4212..77403b6 100644
--- a/src/DAJIN2/core/preprocess/cache_checker.py
+++ b/src/DAJIN2/core/preprocess/cache_checker.py
@@ -2,6 +2,7 @@
 
 import hashlib
 from pathlib import Path
+
 from DAJIN2.utils import io
 
 
diff --git a/src/DAJIN2/core/preprocess/homopolymer_handler.py b/src/DAJIN2/core/preprocess/homopolymer_handler.py
index 1f5ab81..4259b4d 100644
--- a/src/DAJIN2/core/preprocess/homopolymer_handler.py
+++ b/src/DAJIN2/core/preprocess/homopolymer_handler.py
@@ -37,7 +37,7 @@ def extract_sequence_errors_in_homopolymer_loci(
     indels_normalized_control: dict[str, np.array],
     anomal_loci: dict[set],
 ) -> dict[str, set[int]]:
-    sequence_errors_in_homopolymer = dict()
+    sequence_errors_in_homopolymer = {}
     for mut in ["+", "-", "*"]:
         repeat_regions = get_repeat_regions(sequence, anomal_loci[mut])
         if len(repeat_regions) == 0:
diff --git a/src/DAJIN2/core/preprocess/input_formatter.py b/src/DAJIN2/core/preprocess/input_formatter.py
index 1288508..7ac9322 100644
--- a/src/DAJIN2/core/preprocess/input_formatter.py
+++ b/src/DAJIN2/core/preprocess/input_formatter.py
@@ -1,14 +1,12 @@
 from __future__ import annotations
 
 import uuid
-
-from pathlib import Path
-from dataclasses import dataclass
 from collections import defaultdict
-
-from DAJIN2.utils import io, config, fastx_handler
+from dataclasses import dataclass
+from pathlib import Path
 
 from DAJIN2.core import preprocess
+from DAJIN2.utils import config, fastx_handler, io
 
 
 def parse_arguments(arguments: dict) -> tuple:
@@ -64,7 +62,9 @@ def get_genome_coordinates(genome_urls: dict, fasta_alleles: dict, is_cache_geno
         if is_cache_genome:
             genome_coordinates = next(io.read_jsonl(Path(tempdir, "cache", "genome_coordinates.jsonl")))
         else:
-            genome_coordinates = preprocess.fetch_coordinates(genome_coordinates, genome_urls, fasta_alleles["control"])
+            genome_coordinates = preprocess.fetch_coordinates(
+                genome_coordinates, genome_urls, fasta_alleles["control"]
+            )
             genome_coordinates["chrom_size"] = preprocess.fetch_chromosome_size(genome_coordinates, genome_urls)
             io.write_jsonl([genome_coordinates], Path(tempdir, "cache", "genome_coordinates.jsonl"))
 
diff --git a/src/DAJIN2/core/preprocess/insertions_to_fasta.py b/src/DAJIN2/core/preprocess/insertions_to_fasta.py
index dbe95c0..38f572a 100644
--- a/src/DAJIN2/core/preprocess/insertions_to_fasta.py
+++ b/src/DAJIN2/core/preprocess/insertions_to_fasta.py
@@ -1,20 +1,19 @@
 from __future__ import annotations
 
-import uuid
 import random
-from pathlib import Path
+import uuid
+from collections import Counter, defaultdict
 from itertools import groupby
-from collections import defaultdict, Counter
+from pathlib import Path
 from typing import Generator
 
 import numpy as np
 from rapidfuzz import process
 from rapidfuzz.distance import DamerauLevenshtein
-
 from sklearn.cluster import MeanShift
 
 from DAJIN2.core.preprocess.mapping import to_sam
-from DAJIN2.utils import io, config
+from DAJIN2.utils import config, io
 from DAJIN2.utils.cssplits_handler import convert_cssplits_to_cstag
 
 config.set_warnings_ignore()
@@ -65,7 +64,7 @@ def extract_all_insertions(midsv_sample: Generator, mutation_loci: list[set[str]
 def extract_enriched_insertions(
     insertions_sample: dict, insertions_control: dict, coverage_sample: int
 ) -> dict[int, dict[str, int]]:
-    enriched_insertions = dict()
+    enriched_insertions = {}
     threshold_sample = max(5, int(coverage_sample * 0.5 / 100))
     for i in insertions_sample:
         ins_sample: list[str] = insertions_sample[i]
@@ -94,7 +93,7 @@ def extract_enriched_insertions(
             if label in labels_count_sample:
                 del labels_count_sample[label]
 
-        if labels_count_sample == dict():
+        if labels_count_sample == {}:
             continue
 
         # Count the remaining insertion sequences.
@@ -149,7 +148,7 @@ def get_merged_insertion(insertion: dict[str, int], labels: np.ndarray) -> dict[
     insertion_label = [(label, {seq: count}) for label, (seq, count) in zip(labels, insertion.items())]
     insertion_label.sort(key=lambda x: x[0])
 
-    insertion_merged = dict()
+    insertion_merged = {}
     for _, group in groupby(insertion_label, key=lambda x: x[0]):
         group = [g[1] for g in group]
         sequences, counts = set(), 0
@@ -162,7 +161,9 @@ def get_merged_insertion(insertion: dict[str, int], labels: np.ndarray) -> dict[
     return insertion_merged
 
 
-def remove_minor_groups(insertions_merged: dict[tuple[int], dict[tuple[str], int]], coverage: int, percentage: float = 0.5) -> dict[tuple[int], dict[tuple[str], int]]:
+def remove_minor_groups(
+    insertions_merged: dict[tuple[int], dict[tuple[str], int]], coverage: int, percentage: float = 0.5
+) -> dict[tuple[int], dict[tuple[str], int]]:
     """Remove minor groups with less than {percentage} % coverage or less than 5 reads."""
     threshold = max(5, int(coverage * percentage // 100))
     for _, ins in insertions_merged.items():
@@ -184,11 +185,11 @@ def merge_similar_insertions(
 ) -> dict[tuple[int], dict[tuple[str], int]]:
     index_grouped = group_index_by_consecutive_insertions(mutation_loci)
     insertions_grouped = group_insertions(insertions, index_grouped)
-    insertions_merged = dict()
+    insertions_merged = {}
     for idx, insertion in insertions_grouped.items():
         if len(insertion) == 1:
             seq, count = next(iter(insertion.items()))
-            insertions_merged[idx] = {tuple([seq]): count}
+            insertions_merged[idx] = {(seq,): count}
             continue
 
         labels = clustering_insertions(insertion)
@@ -288,7 +289,7 @@ def extract_score_and_sequence(
             scores.append(score)
             sequences.append(",".join(sequence))
 
-    return [(score, sequence) for score, sequence in zip(scores, sequences)]
+    return list(zip(scores, sequences))
 
 
 def filter_minor_label(
@@ -411,7 +412,7 @@ def generate_consensus_insertions(TEMPDIR: Path, SAMPLE_NAME: str, cssplits: lis
 
 def remove_all_n(cons_sequence: dict[int, str]) -> dict[int, str]:
     """Remove all `N` sequences."""
-    cons_sequence_removed = dict()
+    cons_sequence_removed = {}
     for label, seq in cons_sequence.items():
         if all(True if s == "N" else False for s in seq.split(",")):
             continue
@@ -419,9 +420,11 @@ def remove_all_n(cons_sequence: dict[int, str]) -> dict[int, str]:
     return cons_sequence_removed
 
 
-def call_consensus_of_insertion(TEMPDIR: Path, SAMPLE_NAME: str, insertion_sequences_subset: list[dict]) -> dict[int, str]:
+def call_consensus_of_insertion(
+    TEMPDIR: Path, SAMPLE_NAME: str, insertion_sequences_subset: list[dict]
+) -> dict[int, str]:
     """Generate consensus cssplits."""
-    consensus_insertion_cssplits = dict()
+    consensus_insertion_cssplits = {}
     insertion_sequences_subset.sort(key=lambda x: x["LABEL"])
     for label, group in groupby(insertion_sequences_subset, key=lambda x: x["LABEL"]):
         cssplits = [cs["CSSPLIT"].split(",") for cs in group]
@@ -458,7 +461,7 @@ def extract_index_of_insertions(
 def generate_cstag(
     consensus_of_insertions: dict[str, str], index_of_insertions: list[int], sequence: str
 ) -> dict[str, str]:
-    cstag_insertions = dict()
+    cstag_insertions = {}
     for label, cons_seq in consensus_of_insertions.items():
         cons_seq = cons_seq.split(",")
         list_sequence = list(sequence)
@@ -472,7 +475,7 @@ def generate_cstag(
 
 
 def generate_fasta(cstag_insertions: dict[str, str]) -> dict[str, str]:
-    fasta_insertions = dict()
+    fasta_insertions = {}
     for label, cs_tag in cstag_insertions.items():
         fasta_insertions[label] = cstag.to_sequence(cs_tag)
 
@@ -497,7 +500,9 @@ def extract_unique_insertions(fasta_insertions: dict[str, str], FASTA_ALLELES: d
     for key, seq in fasta_insertions_unique.items():
         if key in to_delete:
             continue
-        _, distances, _ = zip(*process.extract_iter(seq, fasta_insertions_unique.values(), scorer=DamerauLevenshtein.distance))
+        _, distances, _ = zip(
+            *process.extract_iter(seq, fasta_insertions_unique.values(), scorer=DamerauLevenshtein.distance)
+        )
         similar_index = {i if d < 10 else None for i, d in enumerate(distances) if i != key}
         to_delete |= similar_index
 
@@ -538,6 +543,7 @@ def save_cstag(TEMPDIR: Path | str, SAMPLE_NAME: str, cstag_insertions: dict[str
     for header, cs_tag in cstag_insertions.items():
         Path(TEMPDIR, SAMPLE_NAME, "cstag", f"{header}.txt").write_text(cs_tag + "\n")
 
+
 ###########################################################
 # main
 ###########################################################
@@ -554,7 +560,7 @@ def generate_insertion_fasta(TEMPDIR, SAMPLE_NAME, CONTROL_NAME, FASTA_ALLELES)
     MUTATION_LOCI = io.load_pickle(Path(TEMPDIR, SAMPLE_NAME, "mutation_loci", "control", "mutation_loci.pickle"))
 
     insertions = extract_insertions(PATH_SAMPLE, PATH_CONTROL, MUTATION_LOCI)
-    if insertions == dict():
+    if insertions == {}:
         """If there is no insertion, return None"""
         return None
 
@@ -578,7 +584,7 @@ def generate_insertion_fasta(TEMPDIR, SAMPLE_NAME, CONTROL_NAME, FASTA_ALLELES)
         [seq for _, seq in insertion_scores_sequences_filtered], labels_filtered, num=1000
     )
     consensus_of_insertions = call_consensus_of_insertion(TEMPDIR, SAMPLE_NAME, insertion_sequences_subset)
-    if consensus_of_insertions == dict():
+    if consensus_of_insertions == {}:
         """
         If there is no insertion sequence, return None
         It is possible when all insertion sequence annotated as `N` that is filtered out
@@ -592,7 +598,7 @@ def generate_insertion_fasta(TEMPDIR, SAMPLE_NAME, CONTROL_NAME, FASTA_ALLELES)
     fasta_insertions = generate_fasta(cstag_insertions)
     fasta_insertions_unique = extract_unique_insertions(fasta_insertions, FASTA_ALLELES)
 
-    if fasta_insertions_unique == dict():
+    if fasta_insertions_unique == {}:
         remove_temporal_files(TEMPDIR, SAMPLE_NAME)
         return None
 
diff --git a/src/DAJIN2/core/preprocess/inversions_to_fasta.py b/src/DAJIN2/core/preprocess/inversions_to_fasta.py
new file mode 100644
index 0000000..1337c77
--- /dev/null
+++ b/src/DAJIN2/core/preprocess/inversions_to_fasta.py
@@ -0,0 +1,615 @@
+from __future__ import annotations
+
+import random
+import uuid
+from collections import Counter, defaultdict
+from itertools import groupby
+from pathlib import Path
+from typing import Generator
+
+import numpy as np
+from rapidfuzz import process
+from rapidfuzz.distance import DamerauLevenshtein
+from sklearn.cluster import MeanShift
+
+from DAJIN2.core.preprocess.mapping import to_sam
+from DAJIN2.utils import config, io
+from DAJIN2.utils.cssplits_handler import convert_cssplits_to_cstag
+
+config.set_warnings_ignore()
+
+import cstag
+
+
+def remove_non_alphabets(cssplit: str) -> str:
+    """Convert a cssplit to a plain DNA sequence."""
+    return "".join([char for char in cssplit if char.isalpha()])
+
+
+###########################################################
+# Cluster insertion alleles
+###########################################################
+
+
+def clustering_insertions(cssplits_insertion: list[str]) -> list[int]:
+    seq_all = [remove_non_alphabets(seq) for seq in cssplits_insertion]
+    query = seq_all[0]
+    _, distances, _ = zip(*process.extract_iter(query, seq_all, scorer=DamerauLevenshtein.normalized_distance))
+
+    insertion_lengths = [[len(c) for c in cs.split(",")] for cs in cssplits_insertion]
+
+    scores = [s + [d] for s, d in zip(insertion_lengths, distances)]
+
+    return MeanShift(bin_seeding=True).fit_predict(np.array(scores)).tolist()
+
+
+###########################################################
+# Detect insertion sequences
+###########################################################
+
+
+def extract_all_inversions(midsv_sample: Generator[list[str]]) -> list[list[int, int, str]]:
+    """Extract sequence of inversions and their start and end indices."""
+    inversions_index_sequence = []
+    for m_sample in midsv_sample:
+        cssplits = m_sample["CSSPLIT"].split(",")
+        min_idx = float("inf")
+        max_idx = -1
+        cs_inversion = []
+        for i, cs in enumerate(cssplits):
+            if cs[-1].islower():
+                min_idx = min(min_idx, i)
+                max_idx = max(max_idx, i)
+                cs_inversion.append(cs)
+        if cs_inversion:
+            inversions_index_sequence.append([min_idx, max_idx, remove_non_alphabets(",".join(cs_inversion))])
+
+    return inversions_index_sequence
+
+
+def cluster_index_of_inversions(
+    inversions_index_sequence: list[list[int, int, str]], coverage: int
+) -> list[tuple[int]]:
+    indexes = [(start, end) for start, end, _ in inversions_index_sequence]
+    min_cluster_size = max(5, int(coverage * 0.5 / 100))
+
+    return MeanShift(bin_seeding=True, min_bin_freq=min_cluster_size).fit_predict(np.array(indexes)).tolist()
+
+
+def extract_enriched_inversions(
+    insertions_sample: dict, insertions_control: dict, coverage_sample: int
+) -> dict[int, dict[str, int]]:
+    enriched_insertions = {}
+    threshold_sample = max(5, int(coverage_sample * 0.5 / 100))
+    for i in insertions_sample:
+        ins_sample: list[str] = insertions_sample[i]
+        ins_control: list[str] = insertions_control.get(i, [])
+
+        labels_all = clustering_insertions(ins_sample + ins_control)
+        labels_sample = labels_all[: len(ins_sample)]
+        labels_control = labels_all[len(ins_sample) : len(ins_sample) + len(ins_control)]
+
+        labels_count_sample = dict(Counter(labels_sample))
+        labels_count_control = dict(Counter(labels_control))
+
+        to_delete = set()
+        # To remove labels containing a high proportion of controls (5% or more, or 5 or more reads).
+        threshold_control = max(5, int(len(labels_control) * 0.05))
+        for label, count_control in labels_count_control.items():
+            if count_control > threshold_control:
+                to_delete.add(label)
+
+        # To delete labels with sample coverage below 0.5% or fewer than 5 reads.
+        for label, count_sample in labels_count_sample.items():
+            if count_sample < threshold_sample:
+                to_delete.add(label)
+
+        for label in to_delete:
+            if label in labels_count_sample:
+                del labels_count_sample[label]
+
+        if labels_count_sample == {}:
+            continue
+
+        # Count the remaining insertion sequences.
+        enriched_insertions[i] = dict(
+            Counter([ins for ins, label in zip(ins_sample, labels_sample) if label in labels_count_sample])
+        )
+    return enriched_insertions
+
+
+def extract_inversions(
+    path_sample: Path, path_control: Path, mutation_loci: list[set[str]]
+) -> dict[int, dict[str, int]]:
+    inversions_sample = extract_all_inversions(io.read_jsonl(path_sample))
+    inversions_control = extract_all_inversions(io.read_jsonl(path_control))
+    coverage_sample = io.count_newlines(path_sample)
+    # coverage_control = io.count_newlines(path_control)
+
+    return extract_enriched_inversions(inversions_sample, inversions_control, coverage_sample)
+
+
+###########################################################
+# Merge similar insertion sequences
+###########################################################
+
+
+def group_index_by_consecutive_insertions(mutation_loci: list[set[str]]) -> list[tuple[int]]:
+    """Groups indices of consecutive locations containing "+" mutations in the input list."""
+    index = sorted(i for i, m in enumerate(mutation_loci) if "+" in m)
+    index_grouped = []
+    for _, group in groupby(enumerate(index), lambda i_x: i_x[0] - i_x[1]):
+        items = [value for _, value in group]
+        index_grouped.append(tuple(items))
+    return index_grouped
+
+
+def group_insertions(
+    insertions: dict[int, dict[str, int]], index_grouped: list[tuple[int]]
+) -> dict[tuple[int], dict[str, int]]:
+    """
+    Insertion bases in consecutive insertion base positions are grouped together in mutation_loci,
+    as the insertion site may shift by a few bases.
+    """
+    insertions_grouped = defaultdict(dict)
+    for idx in index_grouped:
+        for i in idx:
+            if i not in insertions:
+                continue
+            insertions_grouped[idx].update(insertions[i])
+    return dict(insertions_grouped)
+
+
+def get_merged_insertion(insertion: dict[str, int], labels: np.ndarray) -> dict[tuple[str], int]:
+    insertion_label = [(label, {seq: count}) for label, (seq, count) in zip(labels, insertion.items())]
+    insertion_label.sort(key=lambda x: x[0])
+
+    insertion_merged = {}
+    for _, group in groupby(insertion_label, key=lambda x: x[0]):
+        group = [g[1] for g in group]
+        sequences, counts = set(), 0
+        for seq_count in group:
+            seq, count = next(iter(seq_count.items()))
+            sequences.add(seq)
+            counts += count
+        insertion_merged[tuple(sorted(sequences))] = counts
+
+    return insertion_merged
+
+
+def remove_minor_groups(
+    insertions_merged: dict[tuple[int], dict[tuple[str], int]], coverage: int, threshold: float = 0.5
+) -> dict[tuple[int], dict[tuple[str], int]]:
+    for _, ins in insertions_merged.items():
+        # Create a list of elements to delete
+        to_delete = []
+        for seq, count in ins.items():
+            if count < coverage * threshold // 100:
+                to_delete.append(seq)
+
+        # Delete the collected elements
+        for seq in to_delete:
+            del ins[seq]
+
+    return insertions_merged
+
+
+def merge_similar_insertions(
+    insertions: dict[int, dict[str, int]], mutation_loci: list[set[str]], coverage: int, threshold: float = 0.5
+) -> dict[tuple[int], dict[tuple[str], int]]:
+    index_grouped = group_index_by_consecutive_insertions(mutation_loci)
+    insertions_grouped = group_insertions(insertions, index_grouped)
+    insertions_merged = {}
+    for idx, insertion in insertions_grouped.items():
+        if len(insertion) == 1:
+            seq, count = next(iter(insertion.items()))
+            insertions_merged[idx] = {tuple([seq]): count}
+            continue
+
+        labels = clustering_insertions(insertion)
+        insertions_merged[idx] = get_merged_insertion(insertion, labels)
+
+    return remove_minor_groups(insertions_merged, coverage, threshold)
+
+
+###########################################################
+# Pre- and post-process of clustering insertion alleles
+###########################################################
+
+
+def flatten_keys_to_set(dict_keys: list[tuple[int]]) -> set[int]:
+    flattened_nums = set()
+    for key_tuple in dict_keys:
+        for num in key_tuple:
+            flattened_nums.add(num)
+    return flattened_nums
+
+
+def subset_insertions(
+    insertions_merged: dict[tuple[int], dict[tuple[str], int]], num_subset: int = 100
+) -> dict[tuple[int], dict[tuple[str], int]]:
+    """
+    If the number of seqs exceeds `num_subset`, limit it to only the `num_subset` elements with fiexd random sampling.
+    """
+    random.seed(1)
+    insertions_merged_subset = defaultdict(dict)
+    for idx_grouped, insertions in insertions_merged.items():
+        for cs_insertion, counts in insertions.items():
+            if len(cs_insertion) > num_subset:
+                sequences = random.sample(cs_insertion, num_subset)
+            else:
+                sequences = cs_insertion
+            sequences = tuple(remove_non_alphabets(seq) for seq in sequences)
+            insertions_merged_subset[idx_grouped][sequences] = counts
+
+    return dict(insertions_merged_subset)
+
+
+def extract_score_and_sequence(
+    path_sample, insertions_merged: dict[tuple[int], dict[tuple[str], int]]
+) -> list[tuple[list[int], str]]:
+    scores = []
+    sequences = []
+
+    insertions_merged_subset = subset_insertions(insertions_merged, num_subset=100)
+    set_keys = flatten_keys_to_set(insertions_merged.keys())
+
+    cache_score = defaultdict(int)
+    for m in io.read_jsonl(path_sample):
+        cssplits = m["CSSPLIT"].split(",")
+
+        # Skip if there is no insertion
+        if all(not cssplits[i].startswith("+") for i in set_keys):
+            continue
+
+        score = [0] * len(insertions_merged)
+        sequence = ["N"] * len(insertions_merged)
+        for i, (idx_grouped, insertions) in enumerate(insertions_merged.items()):
+            for idx in idx_grouped:
+                if not cssplits[idx].startswith("+"):
+                    continue
+
+                # Use the cache if there is already a precomputed value
+                if (idx, cssplits[idx]) in cache_score:
+                    score[i] = cache_score[(idx, cssplits[idx])]
+                    sequence[i] = cssplits[idx]
+                    continue
+
+                # Record the insertion sequence and count if the sample and insertion sequences are the same
+                for seqs, count in insertions.items():
+                    if cssplits[idx] in seqs:
+                        score[i] = count
+                        sequence[i] = cssplits[idx]
+                        cache_score[(idx, cssplits[idx])] = count
+                        continue
+
+                flag_break = False
+                for seqs, count in insertions_merged_subset[idx_grouped].items():
+                    for _, distance, _ in process.extract_iter(
+                        remove_non_alphabets(cssplits[idx]), seqs, scorer=DamerauLevenshtein.normalized_distance
+                    ):
+                        # Record the insertion sequence and count if the sample and insertion sequences are similar
+                        if distance < 0.1:
+                            score[i] = count
+                            sequence[i] = cssplits[idx]
+                            flag_break = True
+                            break
+
+                    if flag_break:
+                        cache_score[(idx, cssplits[idx])] = count
+                        break
+
+        if any(score):
+            scores.append(score)
+            sequences.append(",".join(sequence))
+
+    return [(score, sequence) for score, sequence in zip(scores, sequences)]
+
+
+def filter_minor_label(
+    labels: list[int],
+    insertions_scores_sequences: list[tuple[list[int], str]],
+    coverage: int,
+    threshold: float = 0.5,
+) -> tuple[list[int], list[str]]:
+    labels_, counts_ = np.unique(labels, return_counts=True)
+    minor_labels = {label for label, count in zip(labels_, counts_) if count < coverage * threshold // 100}
+    index_minor_labels = {i for i, label in enumerate(labels) if label in minor_labels}
+    labels_filtered = [label for i, label in enumerate(labels) if i not in index_minor_labels]
+    score_seq_filterd = [
+        score_seq for i, score_seq in enumerate(insertions_scores_sequences) if i not in index_minor_labels
+    ]
+    return labels_filtered, score_seq_filterd
+
+
+###########################################################
+# Call consensus
+###########################################################
+
+
+def subset_sequences(sequences, labels, num=1000) -> list[dict]:
+    """
+    Downsampling to the required number of sequences when there is a huge number of inserted sequences
+    """
+    sequences_subset = []
+    tmp_sequences = []
+    random.seed(1)
+    for sequence, label in zip(sequences, labels):
+        tmp_sequences.append({"CSSPLIT": sequence, "LABEL": label})
+    tmp_sequences.sort(key=lambda x: x["LABEL"])
+    for _, group in groupby(tmp_sequences, key=lambda x: x["LABEL"]):
+        group = list(group)
+        if len(group) > num:
+            sequences_subset += random.sample(group, num)
+        else:
+            sequences_subset += group
+    return sequences_subset
+
+
+def get_cstag_position(sam_insertions: list[str]) -> tuple[list[str], list[int]]:
+    cs_tags = []
+    positions = []
+    for sam in sam_insertions:
+        if sam.startswith("@"):
+            continue
+        cs_tags.append(sam.split("\t")[-1])
+        positions.append(int(sam.split("\t")[3]))
+    return cs_tags, positions
+
+
+def mapping_insertion(
+    TEMPDIR: Path, SAMPLE_NAME: str, cs_transposed: list[str], consensus_length: int, uuid4: str
+) -> Generator[str]:
+    # Temporarily cache the reference sequence
+    cs_insertion = next(cs for cs in cs_transposed if cs != "N" and len(cs) == consensus_length)
+    ref_seq = cstag.to_sequence(convert_cssplits_to_cstag([cs_insertion]))
+
+    path_reference = Path(TEMPDIR, SAMPLE_NAME, "fasta", f"reference-{uuid4}.fasta")
+    Path(path_reference).write_text(f">reference_insertion\n{ref_seq}\n")
+
+    # Temporarily cache the query sequences
+    path_query = Path(TEMPDIR, SAMPLE_NAME, "fasta", f"query-{uuid4}.fasta")
+    with Path(path_query).open("w") as file:
+        for i, cs_insertion in enumerate(cs_transposed):
+            que_seq = cstag.to_sequence(convert_cssplits_to_cstag([cs_insertion]))
+            file.write(f">query_insertion_{i}\n{que_seq}\n")
+
+    # Switch the mappy settings to either long-read or short-read, depending on the length of the reference sequence
+    if len(ref_seq) > 500:
+        sam_insertions = to_sam(path_reference, path_query)
+    else:
+        options = {"k": 3, "min_dp_score": 1, "min_chain_score": 1, "min_cnt": 1}
+        sam_insertions = to_sam(path_reference, path_query, preset="sr", options=options)
+
+    return sam_insertions
+
+
+def generate_consensus_insertions(TEMPDIR: Path, SAMPLE_NAME: str, cssplits: list[list[str]]) -> str:
+    consensus_insertion = []
+    cssplits_transposed = (list(cs) for cs in zip(*cssplits))
+    for cs_transposed in cssplits_transposed:
+        if all(True if cs == "N" else False for cs in cs_transposed):
+            consensus_insertion.append("N")
+            continue
+
+        count_N = sum(1 for cs in cs_transposed if cs == "N")
+
+        count_cs = defaultdict(int)
+        for cs in cs_transposed:
+            if cs == "N":
+                continue
+            count_cs[len(cs)] += 1
+
+        consensus_length = max(count_cs, key=count_cs.get)
+        count_consensus = sum(1 for cs in cs_transposed if len(cs) == consensus_length)
+
+        if count_N > count_consensus:
+            consensus_insertion.append("N")
+            continue
+
+        sam_insertions = mapping_insertion(TEMPDIR, SAMPLE_NAME, cs_transposed, consensus_length, str(uuid.uuid4()))
+
+        cs_tags, positions = get_cstag_position(sam_insertions)
+        if cs_tags == []:  # When not mapped at all, append `N`.
+            consensus_insertion.append("N")
+        else:
+            cs_insertion_consensus = cstag.to_sequence(cstag.consensus(cs_tags, positions))
+            cs_last = Counter(cs.split("|")[-1] for cs in cs_transposed if len(cs) == consensus_length).most_common()[
+                0
+            ][0]
+            cssplits_insertion_consensus = "|".join("+" + cs for cs in cs_insertion_consensus[:-1]) + "|" + cs_last
+            consensus_insertion.append(cssplits_insertion_consensus)
+
+    return ",".join(consensus_insertion)
+
+
+def remove_all_n(cons_sequence: dict[int, str]) -> dict[int, str]:
+    """Remove all `N` sequences."""
+    cons_sequence_removed = {}
+    for label, seq in cons_sequence.items():
+        if all(True if s == "N" else False for s in seq.split(",")):
+            continue
+        cons_sequence_removed[label] = seq
+    return cons_sequence_removed
+
+
+def call_consensus_of_insertion(
+    TEMPDIR: Path, SAMPLE_NAME: str, insertion_sequences_subset: list[dict]
+) -> dict[int, str]:
+    """Generate consensus cssplits."""
+    consensus_insertion_cssplits = {}
+    insertion_sequences_subset.sort(key=lambda x: x["LABEL"])
+    for label, group in groupby(insertion_sequences_subset, key=lambda x: x["LABEL"]):
+        cssplits = [cs["CSSPLIT"].split(",") for cs in group]
+        consensus_insertion_cssplits[label] = generate_consensus_insertions(TEMPDIR, SAMPLE_NAME, cssplits)
+
+    return remove_all_n(consensus_insertion_cssplits)
+
+
+###########################################################
+# Generate cstag and FASTA
+###########################################################
+
+
+def extract_index_of_insertions(
+    insertions: dict[int, dict[str, int]], insertions_merged: dict[tuple[int], dict[tuple[str], int]]
+) -> list[int]:
+    """`insertions_merged` contains multiple surrounding indices for a single insertion allele. Among them, select the one index where the insertion allele is most frequent."""
+    index_of_insertions = []
+    max_idx = -1
+    for idx_group in insertions_merged:
+        max_val = -1
+        for idx in idx_group:
+            if idx not in insertions:
+                continue
+            if max_val < sum(insertions[idx].values()):
+                max_val = sum(insertions[idx].values())
+                max_idx = idx
+        if max_idx != -1:
+            index_of_insertions.append(max_idx)
+
+    return index_of_insertions
+
+
+def generate_cstag(
+    consensus_of_insertions: dict[str, str], index_of_insertions: list[int], sequence: str
+) -> dict[str, str]:
+    cstag_insertions = {}
+    for label, cons_seq in consensus_of_insertions.items():
+        cons_seq = cons_seq.split(",")
+        list_sequence = list(sequence)
+        for idx, seq in zip(index_of_insertions, cons_seq):
+            if seq == "N":
+                continue
+            list_sequence[idx] = convert_cssplits_to_cstag([seq]) + "="
+        cstag_insertions[label] = "cs:Z:=" + "".join(list_sequence)
+
+    return cstag_insertions
+
+
+def generate_fasta(cstag_insertions: dict[str, str]) -> dict[str, str]:
+    fasta_insertions = {}
+    for label, cs_tag in cstag_insertions.items():
+        fasta_insertions[label] = cstag.to_sequence(cs_tag)
+
+    return fasta_insertions
+
+
+def extract_unique_insertions(FASTA_ALLELES: dict[str, str], fasta_insertions: dict[str, str]) -> dict[str, str]:
+    """
+    Extract unique insertion alleles if they are dissimilar to the FASTA_ALLELES input by the user.
+    "Unique insertion alleles" are defined as sequences that have a difference of more than 10 bases compared to the sequences in FASTA_ALLELES
+    """
+    to_keep = []
+    for query_key, query_seq in fasta_insertions.items():
+        _, distances, _ = zip(
+            *process.extract_iter(query_seq, FASTA_ALLELES.values(), scorer=DamerauLevenshtein.distance)
+        )
+        if all(d > 10 for d in distances):
+            to_keep.append(query_key)
+
+    return {key: fasta_insertions[key] for key in to_keep if key in fasta_insertions}
+
+
+def update_labels(d: dict, FASTA_ALLELES: dict) -> dict:
+    """
+    Update labels to avoid duplicating user-specified alleles.
+    If the allele 'insertion1' exists in FASTA_ALLELES, increment the digits.
+    (insertion1 -> insertion01 -> insertion001...)
+    """
+    user_defined_alleles = set(FASTA_ALLELES)
+    d_values = list(d.values())
+    len_d = len(d_values)
+    digits_up = 0
+    while True:
+        digits = len(str(len_d)) + digits_up
+        d_updated = {f"insertion{i+1:0{digits}}": seq for i, seq in enumerate(d_values)}
+        if user_defined_alleles.isdisjoint(set(d_updated)):
+            break
+        digits_up += 1
+
+    return d_updated
+
+
+###########################################################
+# Save cstag and fasta
+###########################################################
+
+
+def save_fasta(TEMPDIR: Path | str, SAMPLE_NAME: str, fasta_insertions: dict[str, str]) -> None:
+    for header, seq in fasta_insertions.items():
+        Path(TEMPDIR, SAMPLE_NAME, "fasta", f"{header}.fasta").write_text(f">{header}\n{seq}\n")
+
+
+def save_cstag(TEMPDIR: Path | str, SAMPLE_NAME: str, cstag_insertions: dict[str, str]) -> None:
+    for header, cs_tag in cstag_insertions.items():
+        Path(TEMPDIR, SAMPLE_NAME, "cstag", f"{header}.txt").write_text(cs_tag + "\n")
+
+
+###########################################################
+# main
+###########################################################
+
+
+def remove_temporal_files(TEMPDIR: Path, SAMPLE_NAME: str) -> None:
+    _ = [path.unlink() for path in Path(TEMPDIR, SAMPLE_NAME, "fasta").glob("reference-*.fasta")]
+    _ = [path.unlink() for path in Path(TEMPDIR, SAMPLE_NAME, "fasta").glob("query-*.fasta")]
+
+
+def generate_insertion_fasta(TEMPDIR, SAMPLE_NAME, CONTROL_NAME, FASTA_ALLELES) -> None:
+    PATH_SAMPLE = Path(TEMPDIR, SAMPLE_NAME, "midsv", "control.json")
+    PATH_CONTROL = Path(TEMPDIR, CONTROL_NAME, "midsv", "control.json")
+    MUTATION_LOCI = io.load_pickle(Path(TEMPDIR, SAMPLE_NAME, "mutation_loci", "control.pickle"))
+
+    insertions = extract_inversions(PATH_SAMPLE, PATH_CONTROL, MUTATION_LOCI)
+    if insertions == {}:
+        """If there is no insertion, return None"""
+        return None
+
+    coverage: int = io.count_newlines(PATH_SAMPLE)
+
+    insertions_merged = merge_similar_insertions(insertions, MUTATION_LOCI, coverage, threshold=0.5)
+    if all(True if v == {} else False for v in insertions_merged.values()):
+        """"If all the insertion alleles were minor, return None"""
+        return None
+
+    # Clustering similar insertion alleles
+    insertions_scores_sequences = extract_score_and_sequence(PATH_SAMPLE, insertions_merged)
+    cssplits_insertion = [cssplit for _, cssplit in insertions_scores_sequences]
+    labels = clustering_insertions(cssplits_insertion)
+    labels_filtered, insertion_scores_sequences_filtered = filter_minor_label(
+        labels, insertions_scores_sequences, coverage, threshold=0.5
+    )
+
+    # Consensus calling
+    insertion_sequences_subset = subset_sequences(
+        [seq for _, seq in insertion_scores_sequences_filtered], labels_filtered, num=1000
+    )
+    consensus_of_insertions = call_consensus_of_insertion(TEMPDIR, SAMPLE_NAME, insertion_sequences_subset)
+    if consensus_of_insertions == {}:
+        """
+        If there is no insertion sequence, return None
+        It is possible when all insertion sequence annotated as `N` that is filtered out
+        """
+        remove_temporal_files(TEMPDIR, SAMPLE_NAME)
+        return None
+
+    # Format the insertion sequences as cstag and fasta
+    index_of_insertions = extract_index_of_insertions(insertions, insertions_merged)
+    cstag_insertions = generate_cstag(consensus_of_insertions, index_of_insertions, FASTA_ALLELES["control"])
+    fasta_insertions = generate_fasta(cstag_insertions)
+    fasta_insertions_unique = extract_unique_insertions(FASTA_ALLELES, fasta_insertions)
+
+    if fasta_insertions_unique == {}:
+        remove_temporal_files(TEMPDIR, SAMPLE_NAME)
+        return None
+
+    # Update labels
+    cstag_insertions_update = {key: cstag_insertions[key] for key in fasta_insertions_unique.keys()}
+    cstag_insertions_update = update_labels(cstag_insertions_update, FASTA_ALLELES)
+    fasta_insertions_update = update_labels(fasta_insertions_unique, FASTA_ALLELES)
+
+    save_cstag(TEMPDIR, SAMPLE_NAME, cstag_insertions_update)
+    save_fasta(TEMPDIR, SAMPLE_NAME, fasta_insertions_update)
+
+    remove_temporal_files(TEMPDIR, SAMPLE_NAME)
diff --git a/src/DAJIN2/core/preprocess/knockin_handler.py b/src/DAJIN2/core/preprocess/knockin_handler.py
index d3b86a3..9e8c61f 100644
--- a/src/DAJIN2/core/preprocess/knockin_handler.py
+++ b/src/DAJIN2/core/preprocess/knockin_handler.py
@@ -7,7 +7,6 @@
 
 from DAJIN2.core.preprocess import mapping
 
-
 ###########################################################
 # Consider all mutations are possible in the knockin region
 # For large deletion alleles, the deleted sequence becomes the knock-in region, so all mutations within this region are taken into consideration.
diff --git a/src/DAJIN2/core/preprocess/mapping.py b/src/DAJIN2/core/preprocess/mapping.py
index c74b982..73d6992 100644
--- a/src/DAJIN2/core/preprocess/mapping.py
+++ b/src/DAJIN2/core/preprocess/mapping.py
@@ -1,11 +1,11 @@
 from __future__ import annotations
 
-import cstag
-import mappy
-
 from pathlib import Path
 from typing import Generator
 
+import cstag
+import mappy
+
 from DAJIN2.utils import dna_handler
 
 
@@ -14,7 +14,7 @@ def to_sam(
     path_query_fastx: Path,
     preset: str = "map-ont",
     threads: int = 1,
-    options: dict = {},
+    options: dict = None,
     cslong: bool = True,
 ) -> Generator[str]:
     """Align sequences using mappy and Convert PAF to SAM.
@@ -29,6 +29,9 @@ def to_sam(
     Yields:
         str: SAM formatted alignment.
     """
+    if options is None:
+        options = {}
+
     path_reference_fasta = str(path_reference_fasta)
     path_query_fastx = str(path_query_fastx)
 
@@ -90,8 +93,7 @@ def to_sam(
 
             SAM.append("\t".join(alignment))
 
-    for record in SAM:
-        yield record
+    yield from SAM
 
 
 ########################################################################
@@ -100,8 +102,11 @@ def to_sam(
 
 
 def generate_sam(
-    ARGS, paths_fasta: list[str], mappy_options: dict = {}, is_control: bool = False, is_insertion: bool = False
+    ARGS, paths_fasta: list[str], mappy_options: dict = None, is_control: bool = False, is_insertion: bool = False
 ) -> None:
+    if mappy_options is None:
+        mappy_options = {}
+
     if is_control:
         path_fastq = Path(ARGS.tempdir, ARGS.control_name, "fastq", f"{ARGS.control_name}.fastq.gz")
         name = ARGS.control_name
diff --git a/src/DAJIN2/core/preprocess/midsv_caller.py b/src/DAJIN2/core/preprocess/midsv_caller.py
index e0998ea..c83b2be 100644
--- a/src/DAJIN2/core/preprocess/midsv_caller.py
+++ b/src/DAJIN2/core/preprocess/midsv_caller.py
@@ -1,12 +1,12 @@
 from __future__ import annotations
 
-import midsv
-
+from collections import Counter, defaultdict
 from pathlib import Path
 from typing import Generator
-from collections import Counter, defaultdict
 
-from DAJIN2.utils import io, sam_handler, cssplits_handler
+import midsv
+
+from DAJIN2.utils import cssplits_handler, io, sam_handler
 
 
 def has_inversion_in_splice(CIGAR: str) -> bool:
@@ -88,8 +88,7 @@ def extract_best_alignment_length_from_sam(
 
 
 def transform_to_midsv_format(sam: Generator[list[str]]) -> Generator[list[dict]]:
-    for midsv_sample in midsv.transform(sam, midsv=False, cssplit=True, qscore=False, keep=set(["FLAG"])):
-        yield midsv_sample
+    yield from midsv.transform(sam, midsv=False, cssplit=True, qscore=False, keep={"FLAG"})
 
 
 def replace_internal_n_to_d(midsv_sample: Generator[list[dict]], sequence: str) -> Generator[list[dict]]:
@@ -178,7 +177,7 @@ def convert_consecutive_indels_to_match(cssplit: str) -> str:
         cssplit_reversed[i] = current_cs.split("|")[-1]
 
         # Format deletions
-        for k, insertion in enumerate(insertions, 1):
+        for k, _ in enumerate(insertions, 1):
             cssplit_reversed[i + k] = cssplit_reversed[i + k].replace("-", "=")
 
         i += len(insertions) + 1
@@ -221,7 +220,7 @@ def generate_midsv(ARGS, is_control: bool = False, is_insertion: bool = False) -
             """
             Set the destination for midsv as `barcode02/midsv/insertion1.json` when the sample is barcode02 and the allele is insertion1.
             """
-            path_sam_files = list(Path(ARGS.tempdir, name, "sam", allele).glob(f"*.sam"))
+            path_sam_files = list(Path(ARGS.tempdir, name, "sam", allele).glob("*.sam"))
             path_midsv_output = Path(ARGS.tempdir, name, "midsv", allele, f"{name}.jsonl")
 
         preset_cigar_by_qname = extract_preset_and_cigar_by_qname(path_sam_files)
diff --git a/src/DAJIN2/core/preprocess/mutation_extractor.py b/src/DAJIN2/core/preprocess/mutation_extractor.py
index cb5800d..9ecadc1 100644
--- a/src/DAJIN2/core/preprocess/mutation_extractor.py
+++ b/src/DAJIN2/core/preprocess/mutation_extractor.py
@@ -1,10 +1,10 @@
 from __future__ import annotations
 
-import re
 import bisect
+import re
+from collections import defaultdict
 from pathlib import Path
 from typing import Generator
-from collections import defaultdict
 
 from DAJIN2.utils import config
 
@@ -18,8 +18,8 @@
 import numpy as np
 from sklearn.cluster import MiniBatchKMeans
 
-from DAJIN2.utils import io
 from DAJIN2.core.preprocess.homopolymer_handler import extract_sequence_errors_in_homopolymer_loci
+from DAJIN2.utils import io
 
 
 def count_indels(midsv_sample: Generator[dict], sequence: str) -> dict[str, list[int]]:
@@ -41,7 +41,7 @@ def count_indels(midsv_sample: Generator[dict], sequence: str) -> dict[str, list
 
 
 def normalize_indels(count: dict[str, list[int]]) -> dict[str, np.array]:
-    count_normalized = dict()
+    count_normalized = {}
     match_count = np.array(count["="])
     for mut, indel_count in count.items():
         numerator = np.array(indel_count)
@@ -56,7 +56,7 @@ def minimize_mutation_counts(
     """
     In cases where control has a larger value than sample, adjust the value of sample to match that of control.
     """
-    indels_control_minimized = dict()
+    indels_control_minimized = {}
     for mut in {"+", "-", "*"}:
         indels_control_minimized[mut] = np.minimum(indels_control[mut], indels_sample[mut])
     return indels_control_minimized
@@ -136,7 +136,7 @@ def extract_anomal_loci(
     is_consensus: bool = False,
 ) -> dict[str, set[int]]:
     """Extract outlier loci compareing indel counts between sample and control."""
-    anomal_loci = dict()
+    anomal_loci = {}
     for mut in {"+", "-", "*"}:
         values_sample = indels_normalized_sample[mut]
         values_control = indels_normalized_control[mut]
@@ -171,7 +171,7 @@ def count_elements_within_range(arr, lower_bound, upper_bound):
 
 def merge_index_of_consecutive_indel(mutation_loci: dict[str, set[int]]) -> dict[str, set[int]]:
     """Treat as contiguous indels if there are insertions/deletions within five bases of each other"""
-    mutation_loci_merged = dict()
+    mutation_loci_merged = {}
 
     """Reflect point mutations as they are"""
     mutation_loci_merged["*"] = mutation_loci["*"]
@@ -232,14 +232,14 @@ def summarize_indels(path_midsv: Path, sequence: str) -> tuple:
 
 
 def merge_loci(dissimilar_loci: dict[str, set], anomal_loci: dict[str, set]) -> dict[str, set]:
-    mutation_loci = dict()
+    mutation_loci = {}
     for mut in {"+", "-", "*"}:
         mutation_loci[mut] = dissimilar_loci[mut] | anomal_loci[mut]
     return mutation_loci
 
 
 def add_knockin_loci(candidate_loci: dict[str, set], knockin_loci: set):
-    mutation_loci = dict()
+    mutation_loci = {}
     for mut in {"+", "-", "*"}:
         mutation_loci[mut] = candidate_loci[mut] | knockin_loci
     return mutation_loci
@@ -289,9 +289,11 @@ def extract_mutation_loci(
     path_indels_normalized_sample: Path,
     path_indels_normalized_control: Path,
     path_knockin: Path,
-    thresholds: dict[str, float] = {"*": 0.5, "-": 0.5, "+": 0.5},
+    thresholds: dict[str, float] = None,
     is_consensus: bool = False,
 ) -> list[set[str]]:
+    if thresholds is None:
+        thresholds = {"*": 0.5, "-": 0.5, "+": 0.5}
     indels_normalized_sample = io.load_pickle(path_indels_normalized_sample)
     indels_normalized_control = io.load_pickle(path_indels_normalized_control)
 
diff --git a/src/DAJIN2/core/report/__init__.py b/src/DAJIN2/core/report/__init__.py
index 23934c6..88a13eb 100644
--- a/src/DAJIN2/core/report/__init__.py
+++ b/src/DAJIN2/core/report/__init__.py
@@ -1,3 +1,1 @@
-from DAJIN2.core.report import bam_exporter
-from DAJIN2.core.report import sequence_exporter
-from DAJIN2.core.report import mutation_exporter
+from DAJIN2.core.report import bam_exporter, mutation_exporter, sequence_exporter
diff --git a/src/DAJIN2/core/report/bam_exporter.py b/src/DAJIN2/core/report/bam_exporter.py
index 0d0eae3..7895900 100644
--- a/src/DAJIN2/core/report/bam_exporter.py
+++ b/src/DAJIN2/core/report/bam_exporter.py
@@ -98,7 +98,10 @@ def write_sam_to_bam(sam: list[list[str]], path_sam: str | Path, path_bam: str |
     pysam.index("-@", f"{threads}", str(path_bam))
 
 
-def update_sam(sam: list, GENOME_COODINATES: dict = {}) -> list:
+def update_sam(sam: list, GENOME_COODINATES: dict = None) -> list:
+    if GENOME_COODINATES is None:
+        GENOME_COODINATES = {}
+
     sam_records = sam.copy()
     sam_records = sam_handler.remove_microhomology(sam_records)
     if GENOME_COODINATES["genome"]:
@@ -126,8 +129,8 @@ def export_to_bam(TEMPDIR, NAME, GENOME_COODINATES, THREADS, UUID, RESULT_SAMPLE
     sam_headers = [s for s in sam_updated if s[0].startswith("@")]
     sam_contents = [s for s in sam_updated if not s[0].startswith("@")]
     if is_control:
-        qnames: set[str] = set(list(set(s[0] for s in sam_contents[:10000]))[:100])
-        sam_subset = [s for s in sam_updated if s[0] in qnames]
+        qnames_100reads: set[str] = set(list({s[0] for s in sam_contents[:10000]})[:100])  # subset 100 reads
+        sam_subset = [s for s in sam_updated if s[0] in qnames_100reads]
         path_sam_output = Path(TEMPDIR, "report", "BAM", f"temp_{UUID}_{NAME}_control_cache.sam")
         path_bam_output = Path(TEMPDIR, "cache", ".igvjs", NAME, "control.bam")
         write_sam_to_bam(sam_headers + sam_subset, path_sam_output, path_bam_output, THREADS)
diff --git a/src/DAJIN2/core/report/mutation_exporter.py b/src/DAJIN2/core/report/mutation_exporter.py
index f388243..e799eb5 100644
--- a/src/DAJIN2/core/report/mutation_exporter.py
+++ b/src/DAJIN2/core/report/mutation_exporter.py
@@ -1,12 +1,9 @@
 from __future__ import annotations
 
-
-from pathlib import Path
 from itertools import groupby
-from DAJIN2.utils.cssplits_handler import revcomp_cssplits
+from pathlib import Path
 
-# import midsv
-# from DAJIN2.core import preprocess
+from DAJIN2.utils.cssplits_handler import reallocate_insertion_within_deletion, revcomp_cssplits
 
 ###########################################################
 # group by mutation
@@ -21,29 +18,18 @@ def group_by_mutation(cssplits: list[str]) -> list[list[str]]:
     return [list(group) for _, group in groupby(cssplits, key=lambda x: x[0])]
 
 
-def flatten(lst: list[list]) -> list:
-    results_flattend = []
-    for result in lst:
-        if isinstance(result[0], list):
-            for res in result:
-                results_flattend.append(res)
-        else:
-            results_flattend.append(result)
-    return results_flattend
-
-
 ###########################################################
-# report mutations
+# Report mutations
 ###########################################################
 
 
-def _handle_match(group, genome, start, end, header, chromosome):
+def _handle_match(group, genome, start, end, header, chromosome) -> tuple[list[list[str]], int, int]:
     end += len(group)
     start = end
-    return None, start, end
+    return [None], start, end
 
 
-def _handle_substitution(group, genome, start, end, header, chromosome):
+def _handle_substitution(group, genome, start, end, header, chromosome) -> tuple[list[list[str]], int, int]:
     result = []
     for g in group:
         ref = g[1]
@@ -51,22 +37,20 @@ def _handle_substitution(group, genome, start, end, header, chromosome):
         result.append([header, genome, chromosome, start, end, f"substitution: {ref}>{mut}"])
         end += 1
         start = end
-    if len(result) == 1:
-        result = result[0]
     return result, start, end
 
 
-def _handle_deletion(group, genome, start, end, header, chromosome):
+def _handle_deletion(group, genome, start, end, header, chromosome) -> tuple[list[list[str]], int, int]:
     end += len(group) - 1
     size = len(group)
     seq = "".join([g[-1] for g in group])
     result = [header, genome, chromosome, start, end, f"{size}bp deletion: {seq}"]
     end += 1
     start = end
-    return result, start, end
+    return [result], start, end
 
 
-def _handle_insertion(group, genome, start, end, header, chromosome):
+def _handle_insertion(group, genome, start, end, header, chromosome) -> tuple[list[list[str]], int, int]:
     group = group[0]
     size = group.count("|")
     seq_insertion = "".join([g[-1] for g in group.split("|")[:-1]])
@@ -77,33 +61,31 @@ def _handle_insertion(group, genome, start, end, header, chromosome):
         pass
     elif seq_last.startswith("-"):
         result.append([header, genome, chromosome, start, end, f"1bp deletion: {seq_last[-1]}"])
-    if len(result) == 1:
-        result = result[0]
     end += 1
     start = end
     return result, start, end
 
 
-def _handle_inversion(group, genome, start, end, header, chromosome):
+def _handle_inversion(group, genome, start, end, header, chromosome) -> tuple[list[list[str]], int, int]:
     end += len(group) - 1
     size = len(group)
     seq = "".join([g[-1].upper() for g in group])
     result = [header, genome, chromosome, start, end, f"{size}bp inversion: {seq}"]
     end += 1
     start = end
-    return result, start, end
+    return [result], start, end
 
 
-def _handle_unknown(group, genome, start, end, header, chromosome):
+def _handle_unknown(group, genome, start, end, header, chromosome) -> tuple[list[list[str]], int, int]:
     end += len(group) - 1
     size = len(group)
     result = [header, genome, chromosome, start, end, f"{size}bp unknown bases"]
     end += 1
     start = end
-    return result, start, end
+    return [result], start, end
 
 
-def report_mutations(cssplits_grouped, GENOME_COORDINATES, header):
+def report_mutations(cssplits_grouped: list[list[str]], GENOME_COORDINATES, header) -> list[list[str]]:
     genome = GENOME_COORDINATES["genome"]
     chromosome = GENOME_COORDINATES["chrom"]
     start = end = GENOME_COORDINATES["start"]
@@ -122,8 +104,8 @@ def report_mutations(cssplits_grouped, GENOME_COORDINATES, header):
                 result, start, end = handler(group, genome, start, end, header, chromosome)
                 if prefix == "=":
                     continue
-                results.append(result)
-    return flatten(results)
+                results.extend(result)
+    return [list(map(str, r)) for r in results]
 
 
 ###########################################################
@@ -131,16 +113,19 @@ def report_mutations(cssplits_grouped, GENOME_COORDINATES, header):
 ###########################################################
 
 
-def export_to_csv(TEMPDIR: Path | str, SAMPLE_NAME: str, GENOME_COORDINATES: dict, cons_percentage: dict[list]) -> None:
+def export_to_csv(TEMPDIR: Path, SAMPLE_NAME: str, GENOME_COORDINATES: dict, cons_percentage: dict[str, list]) -> None:
     results = [["Allele ID", "Genome", "Chromosome", "Start", "End", "Mutation"]]
     for header, cons in cons_percentage.items():
         cssplits = [max(c, key=c.get) for c in cons]
         if GENOME_COORDINATES.get("strand") == "-":
             cssplits = revcomp_cssplits(cssplits)
+        cssplits = reallocate_insertion_within_deletion(cssplits, bin_size=500, percentage=50)
         cssplits_inversion = annotate_inversion(cssplits)
         cssplits_grouped = group_by_mutation(cssplits_inversion)
         result = report_mutations(cssplits_grouped, GENOME_COORDINATES, header)
         results.extend(result)
+
     results_csv = "\n".join([",".join(map(str, r)) for r in results]) + "\n"
+
     path_output = Path(TEMPDIR, "report", "MUTATION_INFO", f"{SAMPLE_NAME}.csv")
     path_output.write_text(results_csv)
diff --git a/src/DAJIN2/core/report/sequence_exporter.py b/src/DAJIN2/core/report/sequence_exporter.py
index a0518b5..5bffbd8 100644
--- a/src/DAJIN2/core/report/sequence_exporter.py
+++ b/src/DAJIN2/core/report/sequence_exporter.py
@@ -1,11 +1,12 @@
 from __future__ import annotations
 
-import cstag
 import textwrap
 from pathlib import Path
 
-from DAJIN2.utils.cssplits_handler import convert_cssplits_to_cstag, reallocate_insertion_within_deletion
+import cstag
+
 from DAJIN2.core.report.insertion_reflector import reflect_ref_insertion_to_query
+from DAJIN2.utils.cssplits_handler import convert_cssplits_to_cstag, reallocate_insertion_within_deletion
 
 
 def convert_to_fasta(header: str, sequence: str) -> str:
@@ -42,7 +43,7 @@ def export_to_fasta(TEMPDIR: Path, SAMPLE_NAME: str, cons_sequence: dict) -> Non
 
 def parse_fasta(file_path: Path) -> tuple[str, str]:
     """Parses a FASTA file and returns the header and concatenated sequence."""
-    with open(file_path, "r") as f:
+    with open(file_path) as f:
         lines = f.readlines()
 
     header = lines[0].strip().lstrip(">")
diff --git a/src/DAJIN2/gui.py b/src/DAJIN2/gui.py
index d00d0b7..16dad36 100644
--- a/src/DAJIN2/gui.py
+++ b/src/DAJIN2/gui.py
@@ -4,16 +4,14 @@
 import shutil
 import socket
 import webbrowser
-
-import pandas as pd
-
+from contextlib import closing, redirect_stderr
 from pathlib import Path
 from threading import Timer
-from contextlib import closing, redirect_stderr
 
+import pandas as pd
+from flask import Flask, render_template, request
 from waitress import serve
 from werkzeug.utils import secure_filename
-from flask import Flask, render_template, request
 
 from DAJIN2 import main
 from DAJIN2.utils import config
diff --git a/src/DAJIN2/main.py b/src/DAJIN2/main.py
index 8ae4080..fc2b0b5 100644
--- a/src/DAJIN2/main.py
+++ b/src/DAJIN2/main.py
@@ -9,22 +9,22 @@
 os.environ["NUMEXPR_NUM_THREADS"] = "1"
 
 
-import sys
-import shutil
-import logging
 import argparse
-from pathlib import Path
+import logging
+import shutil
+import sys
 from copy import deepcopy
 from itertools import groupby
+from pathlib import Path
 
 from DAJIN2 import gui, view
 from DAJIN2.core import core
-from DAJIN2.utils import io, config, report_generator, input_validator, multiprocess
+from DAJIN2.utils import config, input_validator, io, multiprocess, report_generator
 
 
 def generate_report(name: str, logger: logging.Logger) -> None:
     report_generator.report(name)
-    logger.info(f"\N{party popper} Finished! Open {config.DAJIN_RESULTS_DIR}/{name} to see the report.")
+    logger.info(f"\N{PARTY POPPER} Finished! Open {config.DAJIN_RESULTS_DIR}/{name} to see the report.")
 
 
 ################################################################################
@@ -36,8 +36,8 @@ def execute_single_mode(arguments: dict[str]):
     # Set logging to export log to stderr and file
     path_logfile = config.get_logfile()
     logger = config.set_logging(path_logfile)
-    logger.info(f"\N{runner} Start running DAJIN2 version {config.DAJIN_VERSION}")
-    logger.info(f"\N{Personal Computer} {' '.join(sys.argv)}")
+    logger.info(f"\N{RUNNER} Start running DAJIN2 version {config.DAJIN_VERSION}")
+    logger.info(f"\N{PERSONAL COMPUTER} {' '.join(sys.argv)}")
 
     # Validate input files
     input_validator.validate_files(arguments["sample"], arguments["control"], arguments["allele"])
@@ -71,7 +71,9 @@ def validate_headers_of_batch_file(headers: set[str], filepath: str) -> None:
         raise ValueError(f'{filepath} must contain "sample", "control", "allele" and "name" in the header')
 
     if not headers.issubset(accepted_headers):
-        raise ValueError(f'Accepted header names of {filepath} are "sample", "control", "allele", "name", or "genome".')
+        raise ValueError(
+            f'Accepted header names of {filepath} are "sample", "control", "allele", "name", or "genome".'
+        )
 
 
 def create_argument_dict(args: dict, cache_urls_genome: dict, is_control: bool) -> dict[str, str]:
@@ -103,7 +105,7 @@ def run_DAJIN2(
             contents.append(args)
 
     # Return a list of unique dictionaries
-    contents_unique = [dict(item) for item in set(frozenset(d.items()) for d in contents)]
+    contents_unique = [dict(item) for item in {frozenset(d.items()) for d in contents}]
 
     contents_unique.sort(key=lambda x: x["sample"])
 
@@ -126,7 +128,7 @@ def execute_batch_mode(arguments: dict[str]):
     validate_headers_of_batch_file(headers, path_batchfile)
 
     # Validate contents and fetch genome urls
-    cache_urls_genome = dict()
+    cache_urls_genome = {}
     records.sort(key=lambda x: x["name"])
     for _, groups in groupby(records, key=lambda x: x["name"]):
         for args in groups:
@@ -144,8 +146,8 @@ def execute_batch_mode(arguments: dict[str]):
         config.reset_logging()
         path_logfile = config.get_logfile()
         logger = config.set_logging(path_logfile)
-        logger.info(f"\N{runner} Start running DAJIN2 version {config.DAJIN_VERSION}")
-        logger.info(f"\N{Personal Computer} {' '.join(sys.argv)}")
+        logger.info(f"\N{RUNNER} Start running DAJIN2 version {config.DAJIN_VERSION}")
+        logger.info(f"\N{PERSONAL COMPUTER} {' '.join(sys.argv)}")
 
         # Run DAJIN2
         run_DAJIN2(groups, cache_urls_genome, is_control=True, num_workers=arguments["threads"])
@@ -181,7 +183,7 @@ def execute():
     ###############################################################################
 
     def batchmode(args):
-        arguments = dict()
+        arguments = {}
         arguments["file"] = args.file
         arguments["threads"] = input_validator.update_threads(int(args.threads))
         arguments["debug"] = args.debug
@@ -232,7 +234,7 @@ def viewmode(args):
             raise ValueError("the following arguments are required: -a/--allele")
         if args.name is None:
             raise ValueError("the following arguments are required: -n/--name")
-        arguments = dict()
+        arguments = {}
         arguments["sample"] = args.sample
         arguments["control"] = args.control
         arguments["allele"] = args.allele
diff --git a/src/DAJIN2/utils/config.py b/src/DAJIN2/utils/config.py
index 8553b67..f146939 100644
--- a/src/DAJIN2/utils/config.py
+++ b/src/DAJIN2/utils/config.py
@@ -1,15 +1,14 @@
 from __future__ import annotations
 
-import sys
+import datetime
 import logging
+import sys
 import warnings
-import datetime
-
 from pathlib import Path
 
 from sklearn.exceptions import ConvergenceWarning
 
-DAJIN_VERSION = "0.5.0"
+DAJIN_VERSION = "0.5.1"
 DAJIN_RESULTS_DIR = Path("DAJIN_Results")
 TEMP_ROOT_DIR = Path(DAJIN_RESULTS_DIR, ".tempdir")
 
diff --git a/src/DAJIN2/utils/cssplits_handler.py b/src/DAJIN2/utils/cssplits_handler.py
index 2ebd27f..c7fbdd5 100644
--- a/src/DAJIN2/utils/cssplits_handler.py
+++ b/src/DAJIN2/utils/cssplits_handler.py
@@ -1,5 +1,7 @@
 from __future__ import annotations
+
 import re
+
 import cstag
 import numpy as np
 import ruptures as rpt
@@ -220,10 +222,10 @@ def _extract_break_points_of_large_deletions(
             signal_subset = signal[max(0, approximate_bp - 50) : min(approximate_bp + 50, len(signal))]
             bp = -1
             for i, sig in enumerate(signal_subset):
-                if key == "start" and sig == False:
+                if key == "start" and not sig:
                     bp = i + max(0, approximate_bp - 50) + window_start
                     break
-                if key == "end" and sig == True:
+                if key == "end" and sig:
                     bp = i + max(0, approximate_bp - 50) + window_start - 1
                     break
             if bp == -1:
diff --git a/src/DAJIN2/utils/fastx_handler.py b/src/DAJIN2/utils/fastx_handler.py
index 58455ac..41763be 100644
--- a/src/DAJIN2/utils/fastx_handler.py
+++ b/src/DAJIN2/utils/fastx_handler.py
@@ -1,12 +1,13 @@
 from __future__ import annotations
 
-import re
 import gzip
-from pathlib import Path
-
 import random
+import re
+import uuid
+from pathlib import Path
 
 import mappy
+import pysam
 from DAJIN2.utils.io import sanitize_name
 
 #################################################
@@ -58,26 +59,40 @@ def is_gzip_file(path_file: str | Path) -> bool:
     try:
         with Path(path_file).open("rb") as f:
             return f.read(2) == b"\x1f\x8b"
-    except IOError:
+    except OSError:
         return False
 
 
-def save_fastq_as_gzip(TEMPDIR: Path, path_fastx: list[Path], barcode: str) -> None:
+def save_bam_files_to_single_fastq(TEMPDIR: Path, path_bam_files: list[Path], sample_name: str) -> None:
+    path_bam_files = [str(path) for path in path_bam_files]
+    path_concatenated_bam = Path(TEMPDIR, sample_name, "fastq", f"tmp_{sample_name}_{str(uuid.uuid4())}.bam")
+    pysam.merge("-f", "-o", str(path_concatenated_bam), *path_bam_files)
+
+    path_concatenated_fastq = Path(TEMPDIR, sample_name, "fastq", f"{sample_name}.fastq.gz")
+    _ = pysam.fastq("-0", str(path_concatenated_fastq), str(path_concatenated_bam))
+
+    path_concatenated_bam.unlink()
+
+
+def save_fastx_files_to_single_fastq(TEMPDIR: Path, path_input_files: list[Path], sample_name: str) -> None:
     """Merge gzip and non-gzip files into a single gzip file."""
-    with gzip.open(Path(TEMPDIR, barcode, "fastq", f"{barcode}.fastq.gz"), "wb") as merged_file:
-        for path_file in path_fastx:
+    with gzip.open(Path(TEMPDIR, sample_name, "fastq", f"{sample_name}.fastq.gz"), "wb") as merged_file:
+        for path_file in path_input_files:
             if is_gzip_file(path_file):
                 with gzip.open(path_file, "rb") as f:
                     merged_file.write(f.read())
             else:
-                with open(path_file, "r") as f:
+                with open(path_file) as f:
                     merged_file.write(f.read().encode())
 
 
-def save_concatenated_fastx(TEMPDIR: Path, path_directory: Path, name: str) -> None:
-    fastx_suffix = {".fa", ".fq", ".fasta", ".fastq", ".fa.gz", ".fq.gz", ".fasta.gz", ".fastq.gz"}
-    path_fastx = [path for path in path_directory.iterdir() if extract_extention(path) in fastx_suffix]
-    save_fastq_as_gzip(TEMPDIR, path_fastx, name)
+def save_inputs_as_single_fastx(TEMPDIR: Path, path_directory: Path, sample_name: str) -> None:
+    file_suffix = {".fa", ".fq", ".fasta", ".fastq", ".fa.gz", ".fq.gz", ".fasta.gz", ".fastq.gz", ".bam"}
+    path_input_files = [path for path in path_directory.iterdir() if extract_extention(path) in file_suffix]
+    if all(extract_extention(path) == ".bam" for path in path_directory.iterdir()):
+        save_bam_files_to_single_fastq(TEMPDIR, path_input_files, sample_name)
+    else:
+        save_fastx_files_to_single_fastq(TEMPDIR, path_input_files, sample_name)
 
 
 #################################################
@@ -90,7 +105,7 @@ def read_lines(path_file) -> list[str]:
         with gzip.open(path_file, "rt") as f:
             return [line.strip() for line in f]
     else:
-        with open(path_file, "r") as f:
+        with open(path_file) as f:
             return [line.strip() for line in f]
 
 
diff --git a/src/DAJIN2/utils/input_validator.py b/src/DAJIN2/utils/input_validator.py
index 4992c6e..5957cf4 100644
--- a/src/DAJIN2/utils/input_validator.py
+++ b/src/DAJIN2/utils/input_validator.py
@@ -1,16 +1,15 @@
 from __future__ import annotations
 
+import hashlib
 import os
-import re
 import ssl
-import hashlib
-
+import xml.etree.ElementTree as ET
 from pathlib import Path
 from urllib.error import URLError
 from urllib.request import urlopen
-import xml.etree.ElementTree as ET
 
 import mappy
+import pysam
 
 ########################################################################
 # To accommodate cases where a user might input negative values or
@@ -30,18 +29,26 @@ def update_threads(threads: int) -> int:
 ########################################################################
 
 
-def validate_file_existence(path_file: str):
+def validate_file_existence(path_file: str) -> None:
     if not Path(path_file).exists():
         raise FileNotFoundError(f"{path_file} is not found")
 
 
-def validate_fastq_extension(path_fastq: str):
-    if not re.search(r".fastq$|.fastq.gz$|.fq$|.fq.gz$", path_fastq):
-        raise ValueError(f"{path_fastq} requires extensions either 'fastq', 'fastq.gz', 'fq' or 'fq.gz'")
+def return_file_extension(path_file: Path) -> str:
+    extensions = [".fastq", ".fastq.gz", ".fq", ".fq.gz", ".fasta", ".fasta.gz", ".fa", ".fa.gz", ".bam"]
+    file_suffixes = "".join(path_file.suffixes)  # Combine all suffixes into a single string
+
+    for ext in extensions:
+        if file_suffixes == ext:
+            return ext
+
+    raise ValueError(
+        f"{path_file} requires extensions either .fastq, .fastq.gz, .fq, .fq.gz, .fasta, .fasta.gz, .fa, .fa.gz, or .bam"
+    )
 
 
 # Varidate if the file is in the proper format viewing top 100 lines
-def validate_fastq_content(path_fastq: str):
+def validate_fastq_content(path_fastq: str) -> None:
     try:
         headers, seqs, quals = zip(*[(n, s, q) for i, (n, s, q) in enumerate(mappy.fastx_read(path_fastq)) if i < 100])
         # Remove empty elements
@@ -56,7 +63,7 @@ def validate_fastq_content(path_fastq: str):
         raise ValueError(f"{path_fastq} is not a proper FASTQ format")
 
 
-def validate_fasta_content(path_fasta: str):
+def validate_fasta_content(path_fasta: str, allele_file: bool = False) -> None:
     try:
         headers, seqs = zip(*[(n, s) for n, s, _ in mappy.fastx_read(path_fasta)])
         # Remove empty elements
@@ -71,18 +78,42 @@ def validate_fasta_content(path_fasta: str):
 
     if len(headers) != len(set(headers)):
         raise ValueError(f"{path_fasta} must include unique identifiers")
+
     if len(seqs) != len(set(seqs)):
         raise ValueError(f"{path_fasta} must include unique DNA sequences")
-    if "control" not in headers:
+
+    if allele_file and "control" not in headers:
         raise ValueError(f"One of the headers in the {path_fasta} must be '>control'")
 
 
+def validate_bam_content(path_bam: str) -> None:
+    try:
+        _ = pysam.AlignmentFile(path_bam, "rb")
+    except ValueError:
+        raise ValueError(f"{path_bam} is not a proper BAM format")
+
+
 def validate_files(SAMPLE: str, CONTROL: str, ALLELE: str) -> None:
     for path_directory in [CONTROL, SAMPLE]:
-        _ = [validate_fastq_extension(str(path_fastx)) for path_fastx in Path(path_directory).iterdir()]
-        _ = [validate_fastq_content(str(path_fastx)) for path_fastx in Path(path_directory).iterdir()]
+        extentions = {return_file_extension(path_fastx) for path_fastx in Path(path_directory).iterdir()}
+        if len(extentions) == 1:
+            extention = next(iter(extentions))
+        else:
+            raise ValueError(
+                f"{path_directory} contains multiple extensions. Please check if there are any incorrect files."
+            )
+
+        fasta_extensions = {".fasta", ".fa", ".fasta.gz", ".fa.gz"}
+        fastq_extentions = {".fastq", ".fq", ".fastq.gz", ".fq.gz"}
+        if extention in fasta_extensions:
+            _ = [validate_fasta_content(str(path_fastx)) for path_fastx in Path(path_directory).iterdir()]
+        elif extention in fastq_extentions:
+            _ = [validate_fastq_content(str(path_fastx)) for path_fastx in Path(path_directory).iterdir()]
+        else:
+            _ = [validate_bam_content(str(path_bam)) for path_bam in Path(path_directory).iterdir()]
+
     validate_file_existence(ALLELE)
-    validate_fasta_content(ALLELE)
+    validate_fasta_content(ALLELE, allele_file=True)
 
 
 ########################################################################
diff --git a/src/DAJIN2/utils/io.py b/src/DAJIN2/utils/io.py
index cd52b61..42ce91d 100644
--- a/src/DAJIN2/utils/io.py
+++ b/src/DAJIN2/utils/io.py
@@ -1,19 +1,16 @@
 from __future__ import annotations
 
-import re
 import csv
+import hashlib
 import json
 import pickle
-import hashlib
-
-import wslPath
-
-from pathlib import Path
-
+import re
 from io import BufferedReader
+from pathlib import Path
 from typing import Generator
 
-from openpyxl import load_workbook, Workbook
+import wslPath
+from openpyxl import Workbook, load_workbook
 
 ###########################################################
 # Input/Output
@@ -37,7 +34,7 @@ def save_pickle(data: object, file_path: Path) -> None:
 
 
 def read_jsonl(file_path: str | Path) -> Generator[dict]:
-    with open(file_path, "r") as f:
+    with open(file_path) as f:
         for line in f:
             yield json.loads(line)
 
@@ -100,7 +97,7 @@ def read_xlsx(file_path: str | Path) -> list[dict[str, str]]:
     wb = load_workbook(filename=file_path)
     ws = wb.active
 
-    headers = [cell for cell in next(ws.iter_rows(min_row=1, max_row=1, values_only=True))]
+    headers = list(next(ws.iter_rows(min_row=1, max_row=1, values_only=True)))
 
     records = []
     for row in ws.iter_rows(min_row=2, values_only=True):
@@ -114,8 +111,7 @@ def read_xlsx(file_path: str | Path) -> list[dict[str, str]]:
 
 def read_csv(file_path: str | Path) -> list[dict[str, str]]:
     """Load data from a CSV file."""
-    with open(file_path, "r") as csvfile:
-
+    with open(file_path) as csvfile:
         header = [field.strip() for field in next(csv.reader(csvfile))]
 
         records = []
@@ -125,7 +121,7 @@ def read_csv(file_path: str | Path) -> list[dict[str, str]]:
             if all(element is None for element in row):  # Skip rows with all None values
                 continue
             row_trimmed = [field.strip() for field in row]
-            row_data = {h: v for h, v in zip(header, row_trimmed)}
+            row_data = dict(zip(header, row_trimmed))
             records.append(row_data)
 
         return records
diff --git a/src/DAJIN2/utils/multiprocess.py b/src/DAJIN2/utils/multiprocess.py
index c8864fa..0a2aecf 100644
--- a/src/DAJIN2/utils/multiprocess.py
+++ b/src/DAJIN2/utils/multiprocess.py
@@ -1,15 +1,13 @@
 from __future__ import annotations
 
 import io
-import sys
 import logging
+import sys
 import traceback
-
+from collections.abc import Callable
 from itertools import islice
 from multiprocessing import Process, Queue
-
 from typing import Generator
-from collections.abc import Callable
 
 
 def get_error_message_prefix(arg: dict) -> str:
diff --git a/src/DAJIN2/utils/report_generator.py b/src/DAJIN2/utils/report_generator.py
index fd93dc8..20ea6da 100644
--- a/src/DAJIN2/utils/report_generator.py
+++ b/src/DAJIN2/utils/report_generator.py
@@ -4,7 +4,6 @@
 from pathlib import Path
 
 import plotly.express as px
-
 from DAJIN2.utils import io
 from DAJIN2.utils.config import DAJIN_RESULTS_DIR, TEMP_ROOT_DIR
 
diff --git a/src/DAJIN2/utils/sam_handler.py b/src/DAJIN2/utils/sam_handler.py
index 198f177..5323758 100644
--- a/src/DAJIN2/utils/sam_handler.py
+++ b/src/DAJIN2/utils/sam_handler.py
@@ -1,8 +1,8 @@
 from __future__ import annotations
 
 import re
-
 from itertools import groupby
+
 from DAJIN2.utils.dna_handler import revcomp
 
 
@@ -61,13 +61,15 @@ def remove_overlapped_reads(sam: list[list[str]]) -> list[list[str]]:
 ###########################################################
 
 
-def trim_softclip(CIGAR: str, SEQ: str) -> str:
-    cigar_split = split_cigar(CIGAR)
+def trim_softclip(cigar: str, sequence: str) -> str:
+    if sequence == "*":  # In the case of FASTA format, the qual (= sequence) will be "*".
+        return "*"
+    cigar_split = split_cigar(cigar)
     if cigar_split[0].endswith("S"):
-        SEQ = SEQ[int(cigar_split[0][:-1]) :]
+        sequence = sequence[int(cigar_split[0][:-1]) :]
     if cigar_split[-1].endswith("S"):
-        SEQ = SEQ[: -int(cigar_split[-1][:-1])]
-    return SEQ
+        sequence = sequence[: -int(cigar_split[-1][:-1])]
+    return sequence
 
 
 def check_microhomology(curr_seq_trimmed: str, next_seq_trimmed: str) -> int:
@@ -141,12 +143,12 @@ def trim(alignments: list[list[str]]) -> list[list[str]]:
             alignments[idx + 1][5], num_del = trim_cigar_on_microhomology(next_cigar_splitted, len_microhomology)
             alignments[idx + 1][3] = str(int(next_align[3]) + len_microhomology + num_del)
             alignments[idx + 1][9] = next_seq_trimmed[len_microhomology:]
-            alignments[idx + 1][10] = next_qual_trimmed[len_microhomology:]
+            alignments[idx + 1][10] = next_qual_trimmed[len_microhomology:] if next_qual_trimmed != "*" else "*"
         else:
             alignments[idx][5], num_del = trim_cigar_on_microhomology(next_cigar_splitted, len_microhomology)
             alignments[idx][3] = str(int(alignments[idx][3]) + num_del)
             alignments[idx][9] = next_seq_trimmed[:-len_microhomology]
-            alignments[idx][10] = next_qual_trimmed[:-len_microhomology:]
+            alignments[idx][10] = next_qual_trimmed[:-len_microhomology:] if next_qual_trimmed != "*" else "*"
         idx += 1
 
 
diff --git a/src/DAJIN2/view.py b/src/DAJIN2/view.py
index 3ea2ef6..aa8647a 100644
--- a/src/DAJIN2/view.py
+++ b/src/DAJIN2/view.py
@@ -32,7 +32,7 @@ def execute(name: str):
     env = Environment(loader=FileSystemLoader(path_view, encoding="utf8"))
     template = env.get_template("template_igvjs.html")
     params_genome = {"genome": {"exist": False}}
-    params_reference = dict()
+    params_reference = {}
     path_genome = Path(DIR_IGVJS, "genome_symbol.txt")
     if path_genome.exists():
         path_coordinates = Path(DIR_IGVJS, "genome_coordinates.jsonl")


Test patch:
diff --git a/tests/src/clustering/test_appender.py b/tests/src/clustering/test_appender.py
index 8968ce0..a85bc3e 100644
--- a/tests/src/clustering/test_appender.py
+++ b/tests/src/clustering/test_appender.py
@@ -1,4 +1,4 @@
-from src.DAJIN2.core.clustering.appender import add_readnum, add_percent
+from src.DAJIN2.core.clustering.appender import add_percent, add_readnum
 
 
 def test_add_readnum():
diff --git a/tests/src/clustering/test_kmer_generator.py b/tests/src/clustering/test_kmer_generator.py
index dce464d..caff0d7 100644
--- a/tests/src/clustering/test_kmer_generator.py
+++ b/tests/src/clustering/test_kmer_generator.py
@@ -1,7 +1,7 @@
-import pytest
-
 from typing import Generator
 
+import pytest
+
 from src.DAJIN2.core.clustering.kmer_generator import generate_mutation_kmers
 
 
diff --git a/tests/src/clustering/test_label_merger.py b/tests/src/clustering/test_label_merger.py
index 98f4c7a..1b970d3 100644
--- a/tests/src/clustering/test_label_merger.py
+++ b/tests/src/clustering/test_label_merger.py
@@ -1,6 +1,4 @@
-from DAJIN2.core.clustering.label_merger import merge_mixed_cluster
-from DAJIN2.core.clustering.label_merger import map_clusters_to_previous
-from DAJIN2.core.clustering.label_merger import merge_minor_cluster
+from DAJIN2.core.clustering.label_merger import map_clusters_to_previous, merge_minor_cluster, merge_mixed_cluster
 from DAJIN2.core.clustering.label_updator import relabel_with_consective_order
 
 
diff --git a/tests/src/consensus/test_consensus.py b/tests/src/consensus/test_consensus.py
index 4d2b86c..9f49a8e 100644
--- a/tests/src/consensus/test_consensus.py
+++ b/tests/src/consensus/test_consensus.py
@@ -1,10 +1,9 @@
 from src.DAJIN2.core.consensus.consensus import (
-    replace_sequence_error,
     adjust_to_100_percent,
     call_percentage,
+    replace_sequence_error,
 )
 
-
 ###########################################################
 # replace_sequence
 ###########################################################
diff --git a/tests/src/consensus/test_name_handler.py b/tests/src/consensus/test_name_handler.py
index 81dd808..84c0a56 100644
--- a/tests/src/consensus/test_name_handler.py
+++ b/tests/src/consensus/test_name_handler.py
@@ -1,15 +1,14 @@
-import pytest
-
 from collections import defaultdict
 from typing import NamedTuple
 
+import pytest
 from DAJIN2.core.consensus.name_handler import (
     _detect_sv,
-    _format_allele_label,
     _determine_suffix,
+    _format_allele_label,
+    add_key_by_allele_name,
     call_allele_name,
     update_key_by_allele_name,
-    add_key_by_allele_name,
 )
 
 ###########################################################
@@ -103,7 +102,11 @@ def test_call_allele_name(cons_sequences, cons_percentages, FASTA_ALLELES, thres
 @pytest.mark.parametrize(
     "cons, allele_names, expected_output",
     [
-        ({ConsensusKey("control", 1, 100): "value1", ConsensusKey("control", 2, 100): "value2"}, {1: "name1", 2: "name2"}, {"name1": "value1", "name2": "value2"}),
+        (
+            {ConsensusKey("control", 1, 100): "value1", ConsensusKey("control", 2, 100): "value2"},
+            {1: "name1", 2: "name2"},
+            {"name1": "value1", "name2": "value2"},
+        ),
     ],
 )
 def test_update_key_by_allele_name(cons, allele_names, expected_output):
diff --git a/tests/src/preprocess/test_genome_fetcher.py b/tests/src/preprocess/test_genome_fetcher.py
index 4312ed8..f0aa41b 100644
--- a/tests/src/preprocess/test_genome_fetcher.py
+++ b/tests/src/preprocess/test_genome_fetcher.py
@@ -1,6 +1,7 @@
+from pathlib import Path
+
 import pytest
 
-from pathlib import Path
 from src.DAJIN2.core.preprocess import genome_fetcher
 from src.DAJIN2.utils.input_validator import validate_genome_and_fetch_urls
 
diff --git a/tests/src/preprocess/test_mapping.py b/tests/src/preprocess/test_mapping.py
index 1a3228b..c2bb386 100644
--- a/tests/src/preprocess/test_mapping.py
+++ b/tests/src/preprocess/test_mapping.py
@@ -1,6 +1,7 @@
 from __future__ import annotations
 
 from pathlib import Path
+
 from DAJIN2.core.preprocess import mapping
 
 
diff --git a/tests/src/preprocess/test_midsv_caller.py b/tests/src/preprocess/test_midsv_caller.py
index 0b945e5..56d44fb 100644
--- a/tests/src/preprocess/test_midsv_caller.py
+++ b/tests/src/preprocess/test_midsv_caller.py
@@ -1,14 +1,15 @@
 from __future__ import annotations
 
-import pytest
-
-from DAJIN2.core.preprocess.midsv_caller import has_inversion_in_splice
-from DAJIN2.core.preprocess.midsv_caller import replace_internal_n_to_d
-from DAJIN2.core.preprocess.midsv_caller import convert_flag_to_strand
-from DAJIN2.core.preprocess.midsv_caller import convert_consecutive_indels_to_match
-
 from pathlib import Path
+
+import pytest
 from DAJIN2.core.preprocess.mapping import to_sam
+from DAJIN2.core.preprocess.midsv_caller import (
+    convert_consecutive_indels_to_match,
+    convert_flag_to_strand,
+    has_inversion_in_splice,
+    replace_internal_n_to_d,
+)
 
 
 ###########################################################
diff --git a/tests/src/report/test_bam_exporter.py b/tests/src/report/test_bam_exporter.py
index fb26938..fbbbaf4 100644
--- a/tests/src/report/test_bam_exporter.py
+++ b/tests/src/report/test_bam_exporter.py
@@ -1,8 +1,10 @@
 from __future__ import annotations
 
-import pytest
 from pathlib import Path
 from unittest.mock import patch
+
+import pytest
+
 from src.DAJIN2.core.report import bam_exporter
 
 ###############################################################################
diff --git a/tests/src/report/test_mutation_exporter.py b/tests/src/report/test_mutation_exporter.py
index 93fc7b7..9d13339 100644
--- a/tests/src/report/test_mutation_exporter.py
+++ b/tests/src/report/test_mutation_exporter.py
@@ -120,7 +120,7 @@ def test_report_mutations_substitution():
     cssplits_inversion = report.mutation_exporter.annotate_inversion(cssplits)
     cssplits_grouped = report.mutation_exporter.group_by_mutation(cssplits_inversion)
     test = report.mutation_exporter.report_mutations(cssplits_grouped, GENOME_COODINATES, header)
-    answer = [["test", "mm10", "chr1", 1, 1, "substitution: A>G"]]
+    answer = [["test", "mm10", "chr1", "1", "1", "substitution: A>G"]]
     assert test == answer
 
 
@@ -131,7 +131,10 @@ def test_report_mutations_consecutive_substitution():
     cssplits_inversion = report.mutation_exporter.annotate_inversion(cssplits)
     cssplits_grouped = report.mutation_exporter.group_by_mutation(cssplits_inversion)
     test = report.mutation_exporter.report_mutations(cssplits_grouped, GENOME_COODINATES, header)
-    answer = [["test", "mm10", "chr1", 1, 1, "substitution: A>G"], ["test", "mm10", "chr1", 2, 2, "substitution: C>T"]]
+    answer = [
+        ["test", "mm10", "chr1", "1", "1", "substitution: A>G"],
+        ["test", "mm10", "chr1", "2", "2", "substitution: C>T"],
+    ]
     assert test == answer
 
 
@@ -142,7 +145,7 @@ def test_report_mutations_deletion():
     cssplits_inversion = report.mutation_exporter.annotate_inversion(cssplits)
     cssplits_grouped = report.mutation_exporter.group_by_mutation(cssplits_inversion)
     test = report.mutation_exporter.report_mutations(cssplits_grouped, GENOME_COODINATES, header)
-    answer = [["test", "mm10", "chr1", 1, 2, "2bp deletion: GT"]]
+    answer = [["test", "mm10", "chr1", "1", "2", "2bp deletion: GT"]]
     assert test == answer
 
 
@@ -153,7 +156,7 @@ def test_report_mutations_insertion():
     cssplits_inversion = report.mutation_exporter.annotate_inversion(cssplits)
     cssplits_grouped = report.mutation_exporter.group_by_mutation(cssplits_inversion)
     test = report.mutation_exporter.report_mutations(cssplits_grouped, GENOME_COODINATES, header)
-    answer = [["test", "mm10", "chr1", 1, 1, "2bp insertion: GT"]]
+    answer = [["test", "mm10", "chr1", "1", "1", "2bp insertion: GT"]]
     assert test == answer
 
 
@@ -164,7 +167,7 @@ def test_report_mutations_insertion_with_substitution():
     cssplits_inversion = report.mutation_exporter.annotate_inversion(cssplits)
     cssplits_grouped = report.mutation_exporter.group_by_mutation(cssplits_inversion)
     test = report.mutation_exporter.report_mutations(cssplits_grouped, GENOME_COODINATES, header)
-    answer = [["test", "mm10", "chr1", 1, 1, "2bp insertion: GT"]]
+    answer = [["test", "mm10", "chr1", "1", "1", "2bp insertion: GT"]]
     assert test == answer
 
 
@@ -175,7 +178,7 @@ def test_report_mutations_inversion():
     cssplits_inversion = report.mutation_exporter.annotate_inversion(cssplits)
     cssplits_grouped = report.mutation_exporter.group_by_mutation(cssplits_inversion)
     test = report.mutation_exporter.report_mutations(cssplits_grouped, GENOME_COODINATES, header)
-    answer = [["test", "mm10", "chr1", 1, 2, "2bp inversion: AT"]]
+    answer = [["test", "mm10", "chr1", "1", "2", "2bp inversion: AT"]]
     assert test == answer
 
 
@@ -187,14 +190,14 @@ def test_report_mutations_various():
     cssplits_grouped = report.mutation_exporter.group_by_mutation(cssplits_inversion)
     test = report.mutation_exporter.report_mutations(cssplits_grouped, GENOME_COODINATES, header)
     answer = [
-        ["test", "mm10", "chr1", 1, 1, "substitution: A>G"],
-        ["test", "mm10", "chr1", 2, 2, "2bp insertion: AA"],
-        ["test", "mm10", "chr1", 2, 2, "1bp deletion: A"],
-        ["test", "mm10", "chr1", 3, 3, "1bp unknown bases"],
-        ["test", "mm10", "chr1", 4, 4, "1bp deletion: G"],
-        ["test", "mm10", "chr1", 5, 5, "substitution: C>G"],
-        ["test", "mm10", "chr1", 6, 6, "substitution: A>G"],
-        ["test", "mm10", "chr1", 7, 8, "2bp inversion: AT"],
+        ["test", "mm10", "chr1", "1", "1", "substitution: A>G"],
+        ["test", "mm10", "chr1", "2", "2", "2bp insertion: AA"],
+        ["test", "mm10", "chr1", "2", "2", "1bp deletion: A"],
+        ["test", "mm10", "chr1", "3", "3", "1bp unknown bases"],
+        ["test", "mm10", "chr1", "4", "4", "1bp deletion: G"],
+        ["test", "mm10", "chr1", "5", "5", "substitution: C>G"],
+        ["test", "mm10", "chr1", "6", "6", "substitution: A>G"],
+        ["test", "mm10", "chr1", "7", "8", "2bp inversion: AT"],
     ]
     assert test == answer
 
@@ -206,5 +209,5 @@ def test_report_mutations_genome_coodinates():
     cssplits_inversion = report.mutation_exporter.annotate_inversion(cssplits)
     cssplits_grouped = report.mutation_exporter.group_by_mutation(cssplits_inversion)
     test = report.mutation_exporter.report_mutations(cssplits_grouped, GENOME_COODINATES, header)
-    answer = [["test", "mm10", "chrX", 101, 101, "substitution: A>G"]]
+    answer = [["test", "mm10", "chrX", "101", "101", "substitution: A>G"]]
     assert test == answer
diff --git a/tests/src/utils/test_cssplits_handler.py b/tests/src/utils/test_cssplits_handler.py
index c8e9f3e..cd05d18 100644
--- a/tests/src/utils/test_cssplits_handler.py
+++ b/tests/src/utils/test_cssplits_handler.py
@@ -3,7 +3,6 @@
 import pytest
 from DAJIN2.utils import cssplits_handler
 
-
 ###########################################################
 # find_n_boundaries
 ###########################################################
@@ -124,6 +123,7 @@ def test_get_index_of_large_deletions(cssplits, expected):
 # reallocate_insertion_within_deletion
 ###########################################################
 
+
 @pytest.mark.parametrize(
     "cs, expected",
     [
diff --git a/tests/src/utils/test_fastx_handler.py b/tests/src/utils/test_fastx_handler.py
index 7975d5c..a867e17 100644
--- a/tests/src/utils/test_fastx_handler.py
+++ b/tests/src/utils/test_fastx_handler.py
@@ -1,7 +1,6 @@
 from __future__ import annotations
 
 import pytest
-
 from DAJIN2.utils import fastx_handler
 
 ########################################################################
diff --git a/tests/src/utils/test_input_validator.py b/tests/src/utils/test_input_validator.py
index a670f33..edf789a 100644
--- a/tests/src/utils/test_input_validator.py
+++ b/tests/src/utils/test_input_validator.py
@@ -1,11 +1,12 @@
 from __future__ import annotations
 
+from pathlib import Path
+
 import pytest
 from DAJIN2.utils import input_validator
 
-
 ###############################################################################
-# validate File existance
+# validate File existance and the extentions
 ###############################################################################
 
 
@@ -16,18 +17,19 @@ def test_exists():
     assert str(e.value) == f"{test} is not found"
 
 
+def test_return_file_extension():
+    with pytest.raises(ValueError) as e:
+        test = Path("test.fqq")
+        expected = f"{test} requires extensions either .fastq, .fastq.gz, .fq, .fq.gz, .fasta, .fasta.gz, .fa, .fa.gz, or .bam"
+        input_validator.return_file_extension(test)
+    assert str(e.value) == expected
+
+
 ###############################################################################
 # validate FASTQ
 ###############################################################################
 
 
-def test_fastq_extension():
-    with pytest.raises(ValueError) as e:
-        test = "test.fqq"
-        input_validator.validate_fastq_extension("test.fqq")
-    assert str(e.value) == f"{test} requires extensions either 'fastq', 'fastq.gz', 'fq' or 'fq.gz'"
-
-
 def test_validate_fastq_content_empty():
     with pytest.raises(ValueError):
         fastq_path = "tests/data/utils/validate_inputs/empty.fq"
@@ -74,7 +76,7 @@ def test_fasta_error_duplicated_sequences():
 def test_fasta_error_without_control():
     with pytest.raises(ValueError) as e:
         fasta_path = "tests/data/utils/validate_inputs/no_control.fa"
-        _ = input_validator.validate_fasta_content(fasta_path)
+        _ = input_validator.validate_fasta_content(fasta_path, allele_file=True)
     assert str(e.value) == f"One of the headers in the {fasta_path} must be '>control'"
 
 
diff --git a/tests/src/utils/test_io.py b/tests/src/utils/test_io.py
index 425ea31..0b348bf 100644
--- a/tests/src/utils/test_io.py
+++ b/tests/src/utils/test_io.py
@@ -1,10 +1,11 @@
 from __future__ import annotations
 
-import os
 import json
-import pytest
+import os
 from pathlib import Path
 
+import pytest
+
 from src.DAJIN2.utils import io
 
 
@@ -71,7 +72,7 @@ def test_write_jsonl():
     test_filename = "test_output.json"
     io.write_jsonl(file_path=test_filename, data=data_list)
     # Verify if the file has been written correctly
-    with open(test_filename, "r") as f:
+    with open(test_filename) as f:
         lines = f.readlines()
         assert len(lines) == 2
         loaded_data1 = json.loads(lines[0].strip())
diff --git a/tests/src/utils/test_multiprocess.py b/tests/src/utils/test_multiprocess.py
index 0cccaba..e36915b 100644
--- a/tests/src/utils/test_multiprocess.py
+++ b/tests/src/utils/test_multiprocess.py
@@ -1,12 +1,12 @@
 from __future__ import annotations
 
+import tempfile
+from multiprocessing import Queue
+
 import pytest
 
-import tempfile
 from src.DAJIN2.utils import multiprocess
 
-from multiprocessing import Queue
-
 ###########################################################
 # generate_chunks
 ###########################################################
@@ -125,7 +125,7 @@ def test_multiprocessing_execution():
         multiprocess.run(write_value_to_file, values_to_write, num_workers=3)
 
         # Verify if all values are written correctly
-        with open(file_path, "r") as f:
+        with open(file_path) as f:
             written_values = set(map(int, map(str.strip, f.readlines())))
 
         assert written_values == set(range(1, 11))
diff --git a/tests/src/utils/test_sam_handler.py b/tests/src/utils/test_sam_handler.py
index 57f7b20..f24d567 100644
--- a/tests/src/utils/test_sam_handler.py
+++ b/tests/src/utils/test_sam_handler.py
@@ -1,13 +1,19 @@
 from __future__ import annotations
 
-import midsv
 from pathlib import Path
 
-from DAJIN2.utils.sam_handler import split_cigar
-from DAJIN2.utils.sam_handler import calculate_alignment_length
-from DAJIN2.utils.sam_handler import is_header, is_mapped, is_overlapping, remove_overlapped_reads
-from DAJIN2.utils.sam_handler import remove_microhomology
-from DAJIN2.utils.sam_handler import reverse_flag, revcomp_sam
+import midsv
+from DAJIN2.utils.sam_handler import (
+    calculate_alignment_length,
+    is_header,
+    is_mapped,
+    is_overlapping,
+    remove_microhomology,
+    remove_overlapped_reads,
+    revcomp_sam,
+    reverse_flag,
+    split_cigar,
+)
 
 
 def test_split_cigar():