# Compute QC metrics for each variant and position **URL:** \ **Category:** Hail Query & hailctl\ **Created:** [September 10, 2020, 11:39am UTC](https://discuss.hail.is/t/compute-qc-metrics-for-each-variant-and-position/1641 "2020-09-10T11:39:45Z")\ **Posts on this page:** 6\ **Page:** 1
**Author:** ![SSMK](https://avatars.discourse-cdn.com/v4/letter/s/c0e974/32.png) [@SSMK](https://discuss.hail.is/u/SSMK)\ **Post date:** [September 10, 2020, 11:39am UTC](https://discuss.hail.is/t/compute-qc-metrics-for-each-variant-and-position/1641/1 "2020-09-10T11:39:46Z")
Hello Everyone, I am trying to do something like below. This is how my mental model looks like. ![image](https://canada1.discourse-cdn.com/flex036/uploads/hail/original/1X/a680412313af257a902ae4af2512f157b7c0a6fc.png) The input files have other usual VCF columns such as `INFO` and `FORMAT` fields as well. It’s not shown here. The numbers shown in right hand side output table `variant_qc` metrics are dummy. But I expect the output to be in that format. I had the below approach to solve this a) Read one Input file first b) Group by [variant, chr position] and compute the metrics such as `n_het`, `n_home_ref`,`n_hom_alt` (no need of sample ids as we group by variant and chr position). Should I just use `hl.variant_qc` method here? c) Store the output in a `mt` d) Repeat steps a and b for `file 2` e) combine the output of d with step c. I have few questions a) Is it possible to achieve my expected output using `hail`? How can I do this? b) How can I compute those measures by grouping variants and chr position? What aggregate function should I use to get those measures? I understand I can follow this [link](https://hail.is/docs/0.2/hail.GroupedMatrixTable.html) to use group by rows (variant and chr position)? c) I tried `hl.variant_qc` method on my VCF file\_1 which gave me the below output for count command `(2954429, 1058)` But I also found the same count initially when I read the VCF file `(2954429, 1058)` d) Why isn’t there any reduction in size of the matrix table? I was expecting `variant_qc` output stored in a `mt` to be less in size when compared to the input `mt`. Wouldn’t variant\_qc generate statistics based on each variant and chromosome position? I didn’t see the sample Ids in the variant\_qc, so I was of the understanding that it is grouped by variant and chromosome position. May I know how does that work? ---
**Author:** ![johnc1231](https://yyz2.discourse-cdn.com/flex036/user_avatar/discuss.hail.is/johnc1231/32/286_2.png) [@johnc1231](https://discuss.hail.is/u/johnc1231)\ **Post date:** [September 10, 2020, 1:55pm UTC](https://discuss.hail.is/t/compute-qc-metrics-for-each-variant-and-position/1641/2 "2020-09-10T13:55:56Z")
If these VCFs are both for the same samples, but just represent different chromosomes or something, you can do: ```auto path_1 = "path/to/file1.vcf" path_2 = "path/to/file2.vcf" mt = hl.import_vcf([path_1, path_2]) ``` and it’ll import both at once. You shouldn’t need to use grouping methods, as VCFs don’t look the way you described. A matrix table already has one row per locus/allele pair. Each sample gets its own column. `count` just tells you the number of rows in a `MatrixTable`. `variant_qc` won’t change this. ---
**Author:** ![SSMK](https://avatars.discourse-cdn.com/v4/letter/s/c0e974/32.png) [@SSMK](https://discuss.hail.is/u/SSMK)\ **Post date:** [September 10, 2020, 1:59pm UTC](https://discuss.hail.is/t/compute-qc-metrics-for-each-variant-and-position/1641/3 "2020-09-10T13:59:08Z")
@johnc1231 - So, importing multiple files is like combining all VCF files together? ---
**Author:** ![johnc1231](https://yyz2.discourse-cdn.com/flex036/user_avatar/discuss.hail.is/johnc1231/32/286_2.png) [@johnc1231](https://discuss.hail.is/u/johnc1231)\ **Post date:** [September 10, 2020, 2:01pm UTC](https://discuss.hail.is/t/compute-qc-metrics-for-each-variant-and-position/1641/4 "2020-09-10T14:01:29Z")
Yeah. You should read this if you haven’t: [https://samtools.github.io/hts-specs/VCFv4.2.pdf](https://samtools.github.io/hts-specs/VCFv4.2.pdf) VCF’s don’t look the way you drew them in your Excel table. They have all the samples in one row. ---
**Author:** ![SSMK](https://avatars.discourse-cdn.com/v4/letter/s/c0e974/32.png) [@SSMK](https://discuss.hail.is/u/SSMK)\ **Post date:** [September 11, 2020, 10:21am UTC](https://discuss.hail.is/t/compute-qc-metrics-for-each-variant-and-position/1641/5 "2020-09-11T10:21:39Z")
Hi @johnc1231, Actually what I am trying to achieve is like as shown below I would like to group the records by `variant` and `locus` and compute measures such as `n_het`, `n_hom_ref` etc. So, I tried the below ``` result_mt = (mt.group_rows_by(mt.locus,mt.alleles) .aggregate(n_hom_ref=hl.agg.count_where(mt.GT.is_hom_ref()), n_het = hl.agg.count_where(mt['GT'].is_het()))) ``` However, it produced an incorrect output (due to my code) like below (which isn’t expected). I also tried using `agg.counter` instead of `count.where()` but still incorrect output. ![image](https://canada1.discourse-cdn.com/flex036/uploads/hail/original/1X/706ac5d8cd536f7d01140790b21197ade148c693.png) I expect my output to be like this (at the population level). I don’t want sample level data. | Locus | Alleles | n\_hom\_ref | n\_het | | --- | --- | --- | --- | | 2:1234 | [A,C] | 1012 | 345 | | 3:1234 | [T,TA] | 980 | 654 | | 4:1234 | [C,G] | 321 | 789 | | 5:1234 | [T,G] | 456 | 432 | The above metric numbers shown in table are just dummy and doesn’t carry any meaning May I know whats the how can I transform this the right way? ---
**Author:** ![SSMK](https://avatars.discourse-cdn.com/v4/letter/s/c0e974/32.png) [@SSMK](https://discuss.hail.is/u/SSMK)\ **Post date:** [September 11, 2020, 3:21pm UTC](https://discuss.hail.is/t/compute-qc-metrics-for-each-variant-and-position/1641/6 "2020-09-11T15:21:35Z")
I just found out it can be done using `annotate_rows`. Hence wrote the below code to get the output shown above. If its incorrect or any other better way to write this, please do let me know. I am posting it here for the benefit of others ``` MT2 = mt.annotate_rows(n_hom_ref=hl.agg.count_where(mt.GT.is_hom_ref()), n_het = hl.agg.count_where(mt.GT.is_het()), n_hom_alt = hl.agg.count_where(mt.GT.is_hom_var())) ```