update

.gitignore

@@ -1,3 +1,10 @@
+**.bam
+**.csi
+**.bai 
+**.tbl.gz
+**.pdf
+**.bed.gz
+
 # Byte-compiled / optimized / DLL files
 __pycache__/
 *.py[cod]

Rscripts/UDN-translocation.Rmd

@@ -3,6 +3,7 @@
 SM="UDN318336"
 make_hap_plots(SM, show_imprinted=TRUE)
 ```
+
 ```{r}
 tdf = pdf %>% 
     filter(sample==SM) %>%

isoseq-x-inactivation-affects/add_genes_to_old_data.sh

@@ -0,0 +1,8 @@
+bedtools closest \
+  -a <(zcat ./transcripts.old.data.bed.gz | rg -v "^chrom" | bedtools sort ) \
+  -b ~/projects/fire-figures/data/gencode.v42.annotation_TSS.gff3.gz \
+  -f 0.1 \
+  | grep -vw "\-1" \
+  | sed 's/\tchr.*gene_name=/\t/g' | sed 's/;.*//g' \
+  | bgzip -@ 16 \
+  > ./transcripts.old.data.genes.bed.gz

isoseq-x-inactivation-affects/cmds.sh

@@ -0,0 +1,11 @@
+bedtools makewindows -w 50000 -g ~/assemblies/hg38.analysisSet.fa.fai  > windows.bed
+
+
+bedtools intersect -sorted -c \
+  -a windows.bed \
+  -b H1.bam \
+  | bedtools intersect -sorted -c \
+  -a - \
+  -b H2.bam \
+  > windows.with.counts.bed
+

isoseq-x-inactivation-affects/get-hap-counts-by-gene-from-bam.py

@@ -0,0 +1,24 @@
+import pysam
+import tqdm
+
+
+tags = ["HP", "gn", "in", "tn"]
+print("chrom\tstart\tend\tsample\tHP\tGene\tisoform\ttranscript")
+
+for sample, f in [
+    ("GM12878", "./gm12878.haplotagged.bam"),
+    ("UDN318336", "./mat.pat.old.data.tagged.bam"),
+]:
+    bam = pysam.AlignmentFile(f, threads=16)
+
+    for r in tqdm.tqdm(bam.fetch(until_eof=True)):
+        if r.is_secondary or r.is_supplementary or r.is_unmapped:
+            continue
+        rtn = f"{r.reference_name}\t{r.reference_start}\t{r.reference_end}\t{sample}"
+        for tag in tags:
+            if r.has_tag(tag):
+                res = str(r.get_tag(tag))
+            else:
+                res = "NA"
+            rtn += f"\t{res}"
+        print(rtn)

isoseq-x-inactivation-affects/haplotag.gm12878.sh

@@ -0,0 +1,9 @@
+whatshap haplotag \
+  --ignore-read-groups \
+  ../../phased-bams/GM12878_WGS/GM12878_WGS.haplotagged.vcf.gz \
+  ../../isoseq/GM12878.isoseq.bam \
+  -r ~/assemblies/hg38.analysisSet.fa \
+  | samtools view \
+    -b -@ 16 \
+    -o gm12878.haplotagged.bam --write-index
+

isoseq-x-inactivation-affects/haplotag.old.UDN.data.sh

@@ -0,0 +1,9 @@
+whatshap haplotag \
+  --ignore-read-groups \
+  ../../phased-bams/UDN318336_WGS/UDN318336_WGS.haplotagged.vcf.gz \
+  ./UDN318336.isoseq.bam \
+  -r ~/assemblies/hg38.analysisSet.fa \
+  | samtools view \
+    -b -@ 16 \
+    -o mat.pat.old.data.tagged.bam --write-index
+

profiles/checkpoint/config.yaml

@@ -2,7 +2,7 @@ cluster:
   mkdir -p cluster_logs/{rule} &&
   sbatch
     --account=stergachislab
-    --partition=ckpt
+    --partition=compute
     --nodes=1 
     --ntasks-per-node={threads} 
     --mem={resources.mem_mb}