/
models.py
3958 lines (3338 loc) · 141 KB
/
models.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
import contextlib
import itertools
import os
import shlex
import shutil
import subprocess
import tempfile
import time
from datetime import datetime, timedelta
import json
from collections import defaultdict
import aldjemy
import binning
import wrapt
from itertools import chain
import math
import re
import requests
from bgjobs.plugins import BackgroundJobsPluginPoint
from django.contrib.auth import get_user_model
from django.forms import model_to_dict
from django.utils.html import strip_tags
from sqlalchemy import select, func, and_, delete
import uuid as uuid_object
from postgres_copy import CopyManager
from django.db import models, transaction, connection, utils
from django.db.models import Q
from django.contrib.postgres.fields import ArrayField
from django.contrib.postgres.fields import JSONField
from django.contrib.postgres.indexes import GinIndex
from django.core.exceptions import ValidationError
from django.core.urlresolvers import reverse
from django.dispatch import receiver
from django.conf import settings
from django.db.models.signals import pre_delete
from django.utils import timezone
from projectroles.models import Project
from bgjobs.models import (
BackgroundJob,
JobModelMessageMixin,
LOG_LEVEL_CHOICES,
LOG_LEVEL_ERROR,
LOG_LEVEL_INFO,
)
from projectroles.plugins import get_backend_api
from geneinfo.models import Hgnc, EnsemblToGeneSymbol
from genomicfeatures.models import GeneInterval
#: The SQL Alchemy engine to use
from importer.management.helpers import open_file, tsv_reader
from variants.helpers import SQLALCHEMY_ENGINE
from projectroles.app_settings import AppSettingAPI
app_settings = AppSettingAPI()
#: Django user model.
AUTH_USER_MODEL = getattr(settings, "AUTH_USER_MODEL", "auth.User")
#: The User model to use.
User = get_user_model()
#: Threshold for hom/het ratio for identifying sex.
CHRX_HET_HOM_THRESH = 0.45
#: Threshold for relatedness between parent and children.
THRESH_PARENT = 0.6
#: Threshold for relatedness between siblings.
THRESH_SIBLING = 0.6
#: Pedigree value for male.
PED_MALE = 1
#: Pedigree value for female.
PED_FEMALE = 2
#: Create mapping for chromosome as string to chromosome as integer
CHROMOSOME_STR_TO_CHROMOSOME_INT = {
b: a for a, b in enumerate(list(map(str, range(1, 23))) + ["X", "Y", "MT"], 1)
}
def only_source_name(full_name):
"""Helper function that strips SNAPPY suffixes for samples."""
if full_name.count("-") >= 3:
tokens = full_name.split("-")
return "-".join(tokens[:-3])
else:
return full_name
class CaseAwareProject(Project):
"""A project that is aware of its cases"""
class Meta:
proxy = True
def indices(self, _user=None):
"""Return all registered indices."""
return [p.index for p in self.get_active_smallvariant_cases()]
def pedigree(self, _user=None):
"""Concatenate the pedigrees of project's cases."""
result = []
seen = set()
# Only select cases that have an active variant set.
# TODO Perspectively, we need to distinguish between Small and Structural VariantSets.
for case in self.get_active_smallvariant_cases():
for line in case.pedigree:
if line["patient"] not in seen:
result.append(line)
seen.add((case.name, line["patient"]))
return result
def get_filtered_pedigree_with_samples(self, _user=None):
"""Concatenate the pedigrees of project's cases that have samples."""
result = []
seen = set()
# Only select cases that have an active variant set.
# TODO Perspectively, we need to distinguish between Small and Structural VariantSets.
for case in self.get_active_smallvariant_cases():
for line in case.get_filtered_pedigree_with_samples():
if line["patient"] not in seen:
result.append(line)
seen.add((case.name, line["patient"]))
return result
def get_family_with_filtered_pedigree_with_samples(self, _user=None):
"""Concatenate the pedigrees of project's cases that have samples."""
result = defaultdict(list)
seen = set()
# Only select cases that have an active variant set.
# TODO Perspectively, we need to distinguish between Small and Structural VariantSets.
for case in self.get_active_smallvariant_cases():
for line in case.get_filtered_pedigree_with_samples():
if line["patient"] not in seen:
result[case.name].append(line)
seen.add((case.name, line["patient"]))
return dict(result)
def sample_to_case(self):
"""Compute sample-to-case mapping."""
result = {}
# Only select cases that have an active variant set.
# TODO Perspectively, we need to distinguish between Small and Structural VariantSets.
for case in self.get_active_smallvariant_cases():
for line in case.pedigree:
if line["patient"] not in result:
result[line["patient"]] = case
return result
def chrx_het_hom_ratio(self, sample):
"""Forward to appropriate case"""
case = self.sample_to_case().get(sample)
if not case:
return 0.0
else:
return case.chrx_het_hom_ratio(sample)
def sex_errors(self):
"""Concatenate all contained case's sex errors dicts"""
result = {}
disable_sex_check = app_settings.get_app_setting(
"variants", "disable_pedigree_sex_check", project=self
)
if disable_sex_check:
return result
for case in self.case_set.all():
result.update(case.sex_errors(disable_sex_check))
return result
def sex_errors_to_fix(self):
for case in self.case_set.all():
fix_case = case.sex_errors_variant_stats(lambda x: x)
if fix_case:
return True
return False
def get_case_pks(self):
"""Return PKs for cases."""
return [case.pk for case in self.case_set.all()]
def get_members(self):
"""Return concatenated list of members in ``pedigree``."""
return sorted([x["patient"] for x in self.get_filtered_pedigree_with_samples()])
def get_active_smallvariant_cases(self):
"""Return activate cases."""
return list(self.case_set.filter(smallvariantset__state="active"))
def num_small_vars(self):
"""Return total number of small vars in a project."""
return sum(
case.num_small_vars for case in self.case_set.all() if case.num_small_vars is not None
)
def has_variants_and_variant_sets(self):
return all(case.has_variants_and_variant_set() for case in self.case_set.all())
def casealignmentstats(self):
stats = []
for case in self.case_set.all():
variant_set = case.latest_variant_set
if variant_set:
stats.append(variant_set.casealignmentstats)
return stats
def get_annotation_count(self):
return sum(case.get_annotation_count() for case in self.case_set.all())
def sample_variant_stats(self):
stats = []
for case in self.case_set.all():
variant_set = case.latest_variant_set
if variant_set:
for sample in variant_set.variant_stats.sample_variant_stats.all():
stats.append(sample)
return stats
class SmallVariant(models.Model):
""""Information of a single variant, knows its case."""
#: Genome build
release = models.CharField(max_length=32)
#: Variant coordinates - chromosome
chromosome = models.CharField(max_length=32)
#: Variant coordinates - chromosome (numeric)
chromosome_no = models.IntegerField()
#: Variant coordinates - 1-based start position
start = models.IntegerField()
#: Variant coordinates - end position
end = models.IntegerField()
#: UCSC bin
bin = models.IntegerField()
#: Variant coordinates - reference
reference = models.CharField(max_length=512)
#: Variant coordinates - alternative
alternative = models.CharField(max_length=512)
#: Variant type
var_type = models.CharField(max_length=8)
# #: Link to Case ID.
case_id = models.IntegerField()
#: Link to VariantSet ID.
set_id = models.IntegerField()
#: Miscalleneous information as JSONB.
info = JSONField(default={})
#: Genotype information as JSONB
genotype = JSONField()
#: Number of hom. alt. genotypes
num_hom_alt = models.IntegerField(default=0)
#: Number of hom. ref. genotypes
num_hom_ref = models.IntegerField(default=0)
#: Number of het. genotypes
num_het = models.IntegerField(default=0)
#: Number of hemi alt. genotypes
num_hemi_alt = models.IntegerField(default=0)
#: Number of hemi ref. genotypes
num_hemi_ref = models.IntegerField(default=0)
#: Flag if in clinvar
in_clinvar = models.NullBooleanField()
#: Total ExAC allele frequency
exac_frequency = models.FloatField(null=True)
#: Total ExAC homoyzgous count
exac_homozygous = models.IntegerField(null=True)
#: Total ExAC heterozygous count
exac_heterozygous = models.IntegerField(null=True)
#: Total ExAC hemizgous count
exac_hemizygous = models.IntegerField(null=True)
#: Total thousand genomes frequency count
thousand_genomes_frequency = models.FloatField(null=True)
#: Total thousand genomes homozygous count
thousand_genomes_homozygous = models.IntegerField(null=True)
#: Total thousand genomes heterozygous count
thousand_genomes_heterozygous = models.IntegerField(null=True)
#: Total thousand genomes hemizygous count
thousand_genomes_hemizygous = models.IntegerField(null=True)
#: Total gnomAD exomes frequency
gnomad_exomes_frequency = models.FloatField(null=True)
#: Total gnomAD exomes homozygous count
gnomad_exomes_homozygous = models.IntegerField(null=True)
#: Total gnomAD exomes heterozygous count
gnomad_exomes_heterozygous = models.IntegerField(null=True)
#: Total gnomAD exomes hemizygous count
gnomad_exomes_hemizygous = models.IntegerField(null=True)
#: Total gnomAD genomes frequency
gnomad_genomes_frequency = models.FloatField(null=True)
#: Total gnomAD genomes homozygous count
gnomad_genomes_homozygous = models.IntegerField(null=True)
#: Total gnomAD genomes heterozygous count
gnomad_genomes_heterozygous = models.IntegerField(null=True)
#: Total gnomAD genomes hemizygous count
gnomad_genomes_hemizygous = models.IntegerField(null=True)
#: RefSeq gene ID
refseq_gene_id = models.CharField(max_length=16, null=True)
#: RefSeq transcript ID
refseq_transcript_id = models.CharField(max_length=16, null=True)
#: Flag RefSeq transcript coding
refseq_transcript_coding = models.NullBooleanField()
#: RefSeq HGVS coding sequence
refseq_hgvs_c = models.CharField(max_length=512, null=True)
#: RefSeq HGVS protein sequence
refseq_hgvs_p = models.CharField(max_length=512, null=True)
#: RefSeq variant effect list
refseq_effect = ArrayField(models.CharField(max_length=64), null=True)
#: Distance to next RefSeq exon.
refseq_exon_dist = models.IntegerField(null=True)
#: EnsEMBL gene ID
ensembl_gene_id = models.CharField(max_length=16, null=True)
#: EnsEMBL transcript ID
ensembl_transcript_id = models.CharField(max_length=32, null=True)
#: Flag EnsEMBL transcript coding
ensembl_transcript_coding = models.NullBooleanField()
#: EnsEMBL HGVS coding sequence
ensembl_hgvs_c = models.CharField(max_length=512, null=True)
#: EnsEMBL HGVS protein sequence
ensembl_hgvs_p = models.CharField(max_length=512, null=True)
#: EnsEMBL variant effect list
ensembl_effect = ArrayField(models.CharField(max_length=64, null=True))
#: Distance to next ENSEMBL exon.
ensembl_exon_dist = models.IntegerField(null=True)
#: Allow bulk import
objects = CopyManager()
def get_description(self):
"""Return simple string description of variant"""
return "-".join(
map(str, (self.release, self.chromosome, self.start, self.reference, self.alternative))
)
def __repr__(self):
return "-".join(
map(
str,
(
self.set_id,
self.release,
self.chromosome,
self.start,
self.reference,
self.alternative,
),
)
)
class Meta:
indexes = [
# For query: select all variants for a case.
models.Index(fields=["case_id"]),
# For locating variants by coordiante.
models.Index(fields=["case_id", "chromosome", "bin"]),
# Filter query: the most important thing is to reduce the variants for a case quickly. It's questionable
# how much adding homozygous/frequency really adds here. Adding them back should only done when we
# know that it helps.
GinIndex(fields=["case_id", "refseq_effect"]),
GinIndex(fields=["case_id", "ensembl_effect"]),
models.Index(fields=["case_id", "in_clinvar"]),
# Fast allow-list queries of gene.
models.Index(fields=["case_id", "ensembl_gene_id"]),
models.Index(fields=["case_id", "refseq_gene_id"]),
# For mitochondrial frequency join
models.Index(fields=["case_id", "chromosome_no"]),
]
managed = settings.IS_TESTING
db_table = "variants_smallvariant"
class SmallVariantSummary(models.Model):
"""Summary counts for the small variants.
In the database, this is a materialized view.
"""
#: Genome build
release = models.CharField(max_length=32)
#: Variant coordinates - chromosome
chromosome = models.CharField(max_length=32)
#: Variant coordinates - 1-based start position
start = models.IntegerField()
#: Variant coordinates - end position
end = models.IntegerField()
#: UCSC bin
bin = models.IntegerField()
#: Variant coordinates - reference
reference = models.CharField(max_length=512)
#: Variant coordinates - alternative
alternative = models.CharField(max_length=512)
#: Number of hom. ref. genotypes.
count_hom_ref = models.IntegerField()
#: Number of heterozygous genotypes.
count_het = models.IntegerField()
#: Number of hom. alt. genotypes.
count_hom_alt = models.IntegerField()
#: Number of hemi ref. genotypes.
count_hemi_ref = models.IntegerField()
#: Number of hemi alt. genotypes.
count_hemi_alt = models.IntegerField()
class Meta:
managed = settings.IS_TESTING
db_table = "variants_smallvariantsummary"
def refresh_variants_smallvariantsummary():
"""Refresh the ``SmallVariantSummary`` materialized view."""
with transaction.atomic():
bg_job = BackgroundJob.objects.create(
name='Refreshing small variant summaries (aka "in-house database")',
project=None,
job_type=RefreshSmallVariantSummaryBgJob.spec_name,
user=User.objects.get(username=settings.PROJECTROLES_ADMIN_OWNER),
)
refresh_job = RefreshSmallVariantSummaryBgJob.objects.create(bg_job=bg_job)
with refresh_job.marks():
with connection.cursor() as cursor:
try:
# This will fail if the materialized view is empty.
with transaction.atomic():
cursor.execute(
"REFRESH MATERIALIZED VIEW CONCURRENTLY variants_smallvariantsummary"
)
except utils.NotSupportedError:
with transaction.atomic():
cursor.execute("REFRESH MATERIALIZED VIEW variants_smallvariantsummary")
class CaseManager(models.Manager):
"""Manager for custom table-level Case queries"""
# TODO: properly test searching..
def find(self, search_term, _keywords=None):
"""
Return objects or links matching the query.
:param search_term: Search term (string)
:param _keywords: Optional search keywords as key/value pairs (dict)
:return: Python list of BaseFilesfolderClass objects
"""
objects = super().get_queryset().order_by("name")
objects = objects.filter(
Q(name__iexact=search_term) | Q(search_tokens__icontains=[search_term])
)
return objects
CASE_STATUS_CHOICES = (
("initial", "initial"),
("active", "active"),
("closed-unsolved", "closed as unsolved"),
("closed-uncertain", "closed as uncertain"),
("closed-solved", "closed as solved"),
)
class CoreCase(models.Model):
"""Abstract base class for Case core fields."""
#: Name of the case.
name = models.CharField(max_length=512)
#: Identifier of the index in ``pedigree``.
index = models.CharField(max_length=512)
#: Pedigree information, ``list`` of ``dict`` with the information.
pedigree = JSONField()
#: Note field to summarize the current status
notes = models.TextField(default="", null=True, blank=True)
#: Status field
status = models.CharField(max_length=32, default="initial", choices=CASE_STATUS_CHOICES)
#: Tags field
tags = ArrayField(models.CharField(max_length=32), default=list, null=True, blank=True)
class Meta:
unique_together = (("project", "name"),)
abstract = True
class Case(CoreCase):
"""Stores information about a (germline) case."""
class Meta:
ordering = ("-date_modified",)
indexes = [models.Index(fields=["name"])]
#: DateTime of creation
date_created = models.DateTimeField(auto_now_add=True, help_text="DateTime of creation")
#: DateTime of last modification
date_modified = models.DateTimeField(auto_now=True, help_text="DateTime of last modification")
#: UUID used for identification throughout SODAR.
sodar_uuid = models.UUIDField(
default=uuid_object.uuid4, unique=True, help_text="Case SODAR UUID"
)
#: The number of small variants, ``None`` if no small variants have been imported.
num_small_vars = models.IntegerField(
default=None,
null=True,
blank=True,
verbose_name="Small variants",
help_text="Number of small variants, empty if no small variants have been imported",
)
#: The number of structural variants, ``None`` if no structural variants have been imported.
num_svs = models.IntegerField(
default=None,
null=True,
blank=True,
verbose_name="Structural variants",
help_text="Number of structural variants, empty if no structural variants have been imported",
)
#: The project containing this case.
project = models.ForeignKey(Project, help_text="Project in which this objects belongs")
#: Case manager with custom queries, supporting ``find()`` for the search.
objects = CaseManager()
#: List of additional tokens to search for, for aiding search
search_tokens = ArrayField(
models.CharField(max_length=128, blank=True),
default=list,
db_index=True,
help_text="Search tokens",
)
latest_variant_set = models.ForeignKey(
"SmallVariantSet",
default=None,
blank=True,
null=True,
related_name="case_of_latest_variant_set",
)
latest_structural_variant_set = models.ForeignKey(
"svs.StructuralVariantSet",
default=None,
blank=True,
null=True,
related_name="case_of_latest_structuralvariant_set",
)
def latest_variant_set_id(self):
variant_set = self.latest_variant_set
if variant_set:
return variant_set.id
else:
return -1
def has_variants_and_variant_set(self):
return bool(self.latest_variant_set) and bool(self.num_small_vars)
def latest_structural_variant_set_id(self):
structural_variant_set = self.latest_structural_variant_set
if structural_variant_set:
return structural_variant_set.id
else:
return -1
def has_svs_and_structural_variant_set(self):
return bool(self.latest_structural_variant_set) and bool(self.num_svs)
def days_since_modification(self):
return (timezone.now() - self.date_modified).days
def save(self, *args, **kwargs):
"""Override save() to automatically update ``self.search_tokens``"""
self._update_search_tokens()
super().save(*args, **kwargs)
def _update_search_tokens(self):
"""Force-update ``self.search_tokens``, will enable ``.save()`` call to always save."""
# Get all IDs
self.search_tokens = [self.name] + [x["patient"] for x in self.pedigree if x.get("patient")]
# Remove -N1-DNA1-WES1 etc.
self.search_tokens = [
re.sub(r"-\S+\d+-\S+\d+-[^-]+\d+$", "", x) for x in self.search_tokens
]
# Convert to lower case
self.search_tokens = [x.lower() for x in self.search_tokens]
# Strip non-alphanumeric characters
self.search_tokens = [re.sub(r"[^a-zA-Z0-9]", "", x) for x in self.search_tokens]
def get_sex(self, sample):
"""Return ``int``-value sex for the given ``sample`` in ``pedigree``."""
for line in self.pedigree:
if line["patient"] == sample:
return line["sex"]
def get_absolute_url(self):
"""Return absolute URL for the detail view of this case."""
return reverse(
"variants:case-detail",
kwargs={"project": self.project.sodar_uuid, "case": self.sodar_uuid},
)
def get_filter_url(self):
"""Return absolute URL for the filtration view of this case."""
return reverse(
"variants:case-filter",
kwargs={"project": self.project.sodar_uuid, "case": self.sodar_uuid},
)
def get_background_jobs(self):
"""Return list of ``BackgroundJob`` objects."""
# TODO: need to be more dynamic here?
return BackgroundJob.objects.filter(
Q(variants_exportfilebgjob_related__case=self)
| Q(cadd_submission_bg_job__case=self)
| Q(distiller_submission_bg_job__case=self)
| Q(spanr_submission_bg_job__case=self)
| Q(filter_bg_job__case=self)
)
def get_members(self):
"""Return list of members in ``pedigree``."""
return sorted([x["patient"] for x in self.pedigree])
def get_filtered_pedigree_with_samples(self):
"""Return filtered pedigree lines with members with ``has_gt_entries``."""
# TODO: unit test me
return [x for x in self.pedigree if x["has_gt_entries"]]
def get_family_with_filtered_pedigree_with_samples(self):
"""Concatenate the pedigrees of project's cases that have samples."""
return {self.name: self.get_filtered_pedigree_with_samples()}
def get_members_with_samples(self):
"""Returns names of members that genotype information / samples in imported VCF file."""
# TODO: unit test me
return sorted([x["patient"] for x in self.get_filtered_pedigree_with_samples()])
def get_trio_roles(self):
"""Returns a dict with keys mapping ``index``, ``mother``, ``father`` to pedigree member names if present.
"""
result = {"index": self.index}
for member in self.pedigree:
if member["patient"] == self.index:
if member["father"] != "0":
result["father"] = member["father"]
if member["mother"] != "0":
result["mother"] = member["mother"]
return result
def sex_errors_pedigree(self):
"""Return dict of sample to error messages indicating sex assignment errors that can be derived from the
pedigree information.
Inconsistencies can be determined from father/mother name and sex.
"""
fathers = set([m["father"] for m in self.pedigree])
mothers = set([m["mother"] for m in self.pedigree])
result = {}
for m in self.pedigree:
if m["patient"] in fathers and m["sex"] != PED_MALE:
result[m["patient"]] = ["used as father in pedigree but not male"]
if m["patient"] in mothers and m["sex"] != PED_FEMALE:
result[m["patient"]] = ["used as mother in pedigree not female"]
return result
def chrx_het_hom_ratio(self, sample):
"""Return het./hom. ratio on chrX for ``sample``."""
try:
variant_set = self.latest_variant_set
if not variant_set:
return -1
else:
sample_stats = variant_set.variant_stats.sample_variant_stats.get(
sample_name=sample
)
return sample_stats.chrx_het_hom
except SmallVariantSet.variant_stats.RelatedObjectDoesNotExist:
return -1.0
def sex_errors_variant_stats(self, reporter):
"""Return dict of sample to error messages indicating sex assignment errors that can be derived from
het/hom ratio on chrX.
"""
import_job = ImportVariantsBgJob.objects.filter(case_name=self.name).order_by(
"-date_created"
)
if import_job and not import_job[0].bg_job.status == "done":
return {}
try:
ped_sex = {m["patient"]: m["sex"] for m in self.pedigree}
result = {}
variant_set = self.latest_variant_set
if variant_set:
for sample_stats in variant_set.variant_stats.sample_variant_stats.all():
sample = sample_stats.sample_name
stats_sex = 1 if sample_stats.chrx_het_hom < CHRX_HET_HOM_THRESH else 2
if stats_sex != ped_sex[sample]:
result[sample] = reporter(stats_sex)
return result
except SmallVariantSet.variant_stats.RelatedObjectDoesNotExist:
return {}
def sex_errors_to_fix(self):
return self.sex_errors_variant_stats(lambda x: x)
def sex_errors(self, disable_pedigree_sex_check=None):
"""Returns dict mapping sample to error messages from both pedigree and variant statistics."""
result = {}
if disable_pedigree_sex_check is None:
disable_pedigree_sex_check = app_settings.get_app_setting(
"variants", "disable_pedigree_sex_check", project=self.project
)
if disable_pedigree_sex_check:
return result
for sample, msgs in chain(
self.sex_errors_pedigree().items(),
self.sex_errors_variant_stats(
lambda x: [
"sex from pedigree conflicts with one derived from het/hom ratio on chrX"
]
).items(),
):
result.setdefault(sample, [])
result[sample] += msgs
return result
def rel_errors(self):
"""Returns dict mapping sample to list of relationship errors."""
ped_entries = {m["patient"]: m for m in self.pedigree}
result = {}
import_job = ImportVariantsBgJob.objects.filter(case_name=self.name).order_by(
"-date_created"
)
if import_job and not import_job[0].bg_job.status == "done":
return result
try:
variant_set = self.latest_variant_set
if variant_set:
for rel_stats in variant_set.variant_stats.relatedness.all():
relationship = "other"
if (
ped_entries[rel_stats.sample1]["father"]
== ped_entries[rel_stats.sample2]["father"]
and ped_entries[rel_stats.sample1]["mother"]
== ped_entries[rel_stats.sample2]["mother"]
and ped_entries[rel_stats.sample1]["father"] != "0"
and ped_entries[rel_stats.sample1]["mother"] != "0"
):
relationship = "sibling-sibling"
elif (
ped_entries[rel_stats.sample1]["father"] == rel_stats.sample2
or ped_entries[rel_stats.sample1]["mother"] == rel_stats.sample2
or ped_entries[rel_stats.sample2]["father"] == rel_stats.sample1
or ped_entries[rel_stats.sample2]["mother"] == rel_stats.sample1
):
relationship = "parent-child"
if (
relationship == "sibling-sibling"
and rel_stats.relatedness() < THRESH_SIBLING
) or (
relationship == "parent-child" and rel_stats.relatedness() < THRESH_PARENT
):
for sample in (rel_stats.sample1, rel_stats.sample2):
result.setdefault(sample, []).append(
(
"pedigree shows {} relation for {} and {} but variants show low degree "
"of relatedness"
).format(
relationship,
only_source_name(rel_stats.sample1),
only_source_name(rel_stats.sample2),
)
)
return result
except SmallVariantSet.variant_stats.RelatedObjectDoesNotExist:
return {}
def shortened_notes_text(self, max_chars=50):
"""Shorten ``text`` to ``max_chars`` characters if longer."""
if len(self.notes) > max_chars:
return self.notes[:max_chars] + "..."
else:
return self.notes
def get_annotation_sv_count(self):
variants = set()
for record in self.structural_variant_comments.all():
variants.add((record.chromosome, record.start, record.end, record.sv_type))
for record in self.structural_variant_flags.all():
variants.add((record.chromosome, record.start, record.end, record.sv_type))
return len(variants)
def get_annotation_small_variant_count(self):
variants = set()
for record in self.small_variant_flags.all():
variants.add((record.chromosome, record.start, record.reference, record.alternative))
for record in self.small_variant_comments.all():
variants.add((record.chromosome, record.start, record.reference, record.alternative))
for record in self.acmg_ratings.all():
variants.add((record.chromosome, record.start, record.reference, record.alternative))
return len(variants)
def get_annotation_count(self):
"""Return annotation count."""
return self.get_annotation_sv_count() + self.get_annotation_small_variant_count()
def update_terms(self, terms):
"""Given a dict of individual names to list of terms, ensure that the appropriate ``CasePhenotypeTerms``
records exist.
"""
with transaction.atomic():
self.phenotype_terms.all().delete()
for name, lst in terms.items():
self.phenotype_terms.create(individual=name, terms=lst)
def __str__(self):
"""Return case name as human-readable description."""
return self.name
class CasePhenotypeTerms(models.Model):
"""Phenotype annotation for an individual in a case."""
#: The case that this belongs to.
case = models.ForeignKey(
Case, null=False, related_name="phenotype_terms", help_text="Case for this annotation"
)
#: The name of the individual that this belongs to.
individual = models.CharField(max_length=128, null=False, blank=False, help_text="Individual")
#: A list of HPO, Orphanet, and OMIM terms that the case has been annotated with.
terms = ArrayField(
models.CharField(max_length=128, blank=False),
default=list,
help_text="Phenotype annotation terms with HPO, Orphanet, and OMIM terms",
)
class Meta:
unique_together = (("case", "individual"),)
ordering = ("individual",)
@receiver(pre_delete)
def delete_case_cascaded(sender, instance, **kwargs):
"""Signal handler when attempting to delete a case
Bulk deletes are atomic transactions, including pre/post delete signals.
Comment From their code base in `contrib/contenttypes/fields.py`:
```
if bulk:
# `QuerySet.delete()` creates its own atomic block which
# contains the `pre_delete` and `post_delete` signal handlers.
queryset.delete()
```
"""
if sender == Case:
SmallVariant.objects.filter(case_id=instance.id).delete()
for plugin in BackgroundJobsPluginPoint.get_plugins():
for _, klass in plugin.job_specs.items():
bgjobs = []
if hasattr(klass, "case_id"):
bgjobs = klass.objects.filter(case_id=instance.id)
elif hasattr(klass, "case_name"):
bgjobs = klass.objects.filter(project=instance.project, case_name=instance.name)
for bgjob in bgjobs:
bgjob.bg_job.delete()
class CaseComments(models.Model):
"""Comments associated with a case."""
#: UUID of the job
sodar_uuid = models.UUIDField(
default=uuid_object.uuid4, unique=True, help_text="Case SODAR UUID"
)
#: DateTime of creation
date_created = models.DateTimeField(auto_now_add=True, help_text="DateTime of creation")
#: DateTime of last modification
date_modified = models.DateTimeField(auto_now=True, help_text="DateTime of last modification")
#: Case for the comment
case = models.ForeignKey(
Case, null=False, related_name="case_comments", help_text="Case for this comment"
)
#: User who created the comment
user = models.ForeignKey(
AUTH_USER_MODEL,
on_delete=models.SET_NULL,
blank=True,
null=True,
related_name="case_comments",
help_text="User who created the comment",
)
#: User
comment = models.TextField()
class Meta:
ordering = ["date_created"]
def shortened_text(self, max_chars=50):
"""Shorten ``text`` to ``max_chars`` characters if longer."""
if len(self.comment) > max_chars:
return self.comment[:max_chars] + "..."
else:
return self.comment
def get_absolute_url(self):
return self.case.get_absolute_url() + "#comment-%s" % self.sodar_uuid
class SmallVariantSet(models.Model):
"""A set of small variants associated with a case.
This additional step of redirection is created such that a new set of variants can be imported into the database
without affecting existing variants and without using transactions.
"""
#: DateTime of creation
date_created = models.DateTimeField(auto_now_add=True, help_text="DateTime of creation")
#: DateTime of last modification
date_modified = models.DateTimeField(auto_now=True, help_text="DateTime of last modification")
#: The case that the variants are for.
case = models.ForeignKey(Case, null=False, help_text="The case that this set is for")
#: The state of the variant set.
state = models.CharField(
max_length=16,
choices=(("importing", "importing"), ("active", "active"), ("deleting", "deleting")),
null=False,
blank=False,
)
def cleanup_variant_sets(min_age_hours=12):
"""Cleanup old variant sets."""
variant_sets = list(
SmallVariantSet.objects.filter(
date_created__lte=datetime.now() - timedelta(hours=min_age_hours)
).exclude(state="active")
)
smallvariant_table = aldjemy.core.get_meta().tables["variants_smallvariant"]
for variant_set in variant_sets:
SQLALCHEMY_ENGINE.execute(
smallvariant_table.delete().where(
and_(
smallvariant_table.c.set_id == variant_set.id,
smallvariant_table.c.case_id == variant_set.case.id,
)
)
)
variant_set.delete()
class AnnotationReleaseInfo(models.Model):
"""Model to track the database releases used during annotation of a case.
"""
#: Release of genomebuild
genomebuild = models.CharField(max_length=32, default="GRCh37")
#: Name of imported table
table = models.CharField(max_length=512)
#: Timestamp of import
timestamp = models.DateTimeField(auto_now_add=True, editable=False)
#: Data release
release = models.CharField(max_length=512)
#: Link to case
case = models.ForeignKey(Case)
#: Link to variant set
variant_set = models.ForeignKey(SmallVariantSet)
class Meta:
unique_together = ("genomebuild", "table", "variant_set")
class CaseAlignmentStats(models.Model):
"""Store alignment information about alignment statistics for a case."""
#: Reference to the case.
case = models.ForeignKey(Case, null=False)
#: Reference to the small variant set.
variant_set = models.OneToOneField(SmallVariantSet, null=False)
#: The BAM statistics information. On the top level, there is one entry for each sample. Below this are the keys
#: bamstats, idxstats, min_cov_targets, min_cov_bases. Below bamstats, there is an entry from "samtools stats"
#: metric (as retrieved by ``grep '^SN'`` to the metric value. Below idxstats, the contig name is mapped to
#: a dict with "mapped" and "unmapped" entries. Below min_cov_targets is a dict mapping coverage to the percentage
#: of targets having a minimal coverage of the key value. Below min_cov_bases is a dict mapping coverage to
#: the percentage of *bases* having a minimal coverage of the key value.
bam_stats = JSONField(default={}, null=False)
#: Allow bulk import
objects = CopyManager()
class Meta:
indexes = (models.Index(fields=["case"]),)