# Errors after iterative calls to union\_cols

**URL:** https://discuss.hail.is/t/errors-after-iterative-calls-to-union-cols/790
**Category:** Hail Query & hailctl
**Created:** [November 29, 2018, 12:51am UTC](https://discuss.hail.is/t/errors-after-iterative-calls-to-union-cols/790 "2018-11-29T00:51:17Z")
**Posts on this page:** 4
**Page:** 1

<div class="post-metadata">

### Author: ![nlw](https://avatars.discourse-cdn.com/v4/letter/n/dec6dc/32.png) [@nlw](https://discuss.hail.is/u/nlw)
#### Post date: [November 29, 2018, 12:51am UTC](https://discuss.hail.is/t/errors-after-iterative-calls-to-union-cols/790/1 "2018-11-29T00:51:17Z")

</div>

I want to compute the concordance of the platinum genomes na12878 vs 100 of my own na12878 sequences.  
The `hl.concordance()` function seems like a good choice. However, it requires that you have two MatrixTables with matching sample names (it does not support 1:many comparison).

To get around this requirement, I thought I would:

1. create a MatrixTable with a single sample of PG NA12878 (called pg\_biallelic\_vars\_chr22)
2. iteratively `union_cols()` this to itself 100x (na12878\_pg\_100x)
3. duplicate my\_na12878 sample ids, and key off these in the na12878\_pg\_100x
4. save the 100x file so I could do other things with it (such as concordance)

This is the code I came up with (note that I filtered the PG calls down to biallelic and chr22 only):

```auto
n_samples = 100
na12878_pg_100x = pg_biallelic_vars_chr22
for n in range(1,n_samples):
    na12878_pg_100x = na12878_pg_100x.union_cols(pg_biallelic_vars_chr22)

na12878_pg_100x = na12878_pg_100x.add_col_index()
# we need the keyed sample name to match for the "concordance" calculation, so take it from the real na12878 samples
na12878_pg_100x = na12878_pg_100x.annotate_cols(s2=hl.array(my_na12878.s.take(n_samples))[hl.int(na12878_pg_100x.col_idx)])
na12878_pg_100x = na12878_pg_100x.key_cols_by('s2')
na12878_pg_100x.cols().show(25)
fname = 's3://my-s3-path/samples_one_by_one.mt'
na12878_pg_100x.write(fname, overwrite=True) # <-- this is where it fails
na12878_pg_100x = hl.read_matrix_table(fname)

```

(Note that the code fails with writing no matter if i set `stage_locally` to `True` or `False`.)  
I am able to run the above code with `n_samples=25`, but it fails at 50 or more with the following error:

```auto
FatalError: StackOverflowError: null

Java stack trace:
org.apache.spark.SparkException: Job aborted due to stage failure: Task 0 in stage 81.0 failed 4 times, most recent failure: Lost task 0.3 in stage 81.0 (TID 821, ip-172-19-23-116.ec2.internal, executor 153): java.lang.StackOverflowError
	at org.apache.spark.util.ByteBufferInputStream.read(ByteBufferInputStream.scala:49)
	at java.io.ObjectInputStream$PeekInputStream.read(ObjectInputStream.java:2662)
	at java.io.ObjectInputStream$PeekInputStream.readFully(ObjectInputStream.java:2678)
	at java.io.ObjectInputStream$BlockDataInputStream.readInt(ObjectInputStream.java:3179)
	at java.io.ObjectInputStream.readHandle(ObjectInputStream.java:1683)
	at java.io.ObjectInputStream.readObject0(ObjectInputStream.java:1552)
	at java.io.ObjectInputStream.defaultReadFields(ObjectInputStream.java:2286)
	at java.io.ObjectInputStream.readSerialData(ObjectInputStream.java:2210)
	at java.io.ObjectInputStream.readOrdinaryObject(ObjectInputStream.java:2068)
	at java.io.ObjectInputStream.readObject0(ObjectInputStream.java:1572)
...
	at java.io.ObjectInputStream.readObject0(ObjectInputStream.java:1572)
	at java.io.ObjectInputStream.defaultReadFields(ObjectInputStream.java:2286)
	at java.io.ObjectInputStream.readSerialData(ObjectInputStream.java:2210)
	at java.io.ObjectInputStream.readOrdinaryObject(ObjectInputStream.java:2068)
	at java.io.ObjectInputStream.readObject0(ObjectInputStream.java:1572)
	at java.io.ObjectInputStream.readObject(ObjectInputStream.java:430)
	at scala.collection.immutable.List$SerializationProxy.readObject(List.scala:479)
	at sun.reflect.GeneratedMethodAccessor1.invoke(Unknown Source)
	at sun.reflect.DelegatingMethodAccessorImpl.invoke(DelegatingMethodAccessorImpl.java:43)
	at java.lang.reflect.Method.invoke(Method.java:498)
	at java.io.ObjectStreamClass.invokeReadObject(ObjectStreamClass.java:1170)
	at java.io.ObjectInputStream.readSerialData(ObjectInputStream.java:2177)

Hail version: 0.2-961f76d14f1e
Error summary: StackOverflowError: null

```

So there are two issues to address here:

1. why is this having trouble writing to the file?
2. is there a better way to do the 1:many calculation for the na12878 comparisons.

---

<div class="post-metadata">

### Author: ![tpoterba](https://yyz2.discourse-cdn.com/flex036/user_avatar/discuss.hail.is/tpoterba/32/61_2.png) [@tpoterba](https://discuss.hail.is/u/tpoterba)
#### Post date: [November 29, 2018, 1:43am UTC](https://discuss.hail.is/t/errors-after-iterative-calls-to-union-cols/790/2 "2018-11-29T01:43:56Z")

</div>

1. why is this having trouble writing to the file?

It’s important to note that almost every operation in Hail is lazy – that is, everything upstream is executed when you are executing the `write`.

1. is there a better way to do the 1:many calculation for the na12878 comparisons.

I think the iterative union\_cols is quadratic in the number of samples – see this thread involving the same issue:

> [@Performance of writing matrixtable on 0.2](http://discuss.hail.is/t/performance-of-writing-matrixtable-on-0-2/663/9):
>
> for now you can also get down to O(N \* log2 N) doing a tree union: mts\_ = mts[:] iteration = 0 while (len(mts\_) \> 1): iteration += 1 print(f'iteration {iteration}') tmp = [] for i in range(0, len(mts\_), 2): tmp.append(hl.MatrixTable.union\_cols(mts\_[i:i+2])) mts\_ = tmp[:] [final\_mt] = mts\_

The issue is that Hail doesn’t have great ways right now to load a batch of single-sample VCFs into a single joint MatrixTable. We’re working on a gVCF merging algorithm, which we expect to come online in the next six months.

For the time being, things may work if you use my N \* log(N) solution in the linked comment.

---

<div class="post-metadata">

### Author: ![danking](https://yyz2.discourse-cdn.com/flex036/user_avatar/discuss.hail.is/danking/32/43_2.png) [@danking](https://discuss.hail.is/u/danking)
#### Post date: [December 3, 2018, 4:39pm UTC](https://discuss.hail.is/t/errors-after-iterative-calls-to-union-cols/790/3 "2018-12-03T16:39:06Z")

</div>

Would this work?

```auto
In [4]: mt = hl.balding_nichols_model(1, 1, 100)
2018-12-03 11:33:00 Hail: INFO: balding_nichols_model: generating genotypes for 1 populations, 1 samples, and 100 variants...

In [5]: mt.count_cols()
Out[5]: 1

In [15]: mt2 = mt.annotate_cols(s = hl.range(0, 100).map(lambda i: hl.str(mt.sample_idx + i)))

In [16]: mt2 = mt2.explode_cols('s').key_cols_by('s')

In [17]: mt2.count_cols()
2018-12-03 11:37:52 Hail: INFO: Coerced sorted dataset
2018-12-03 11:37:53 Hail: INFO: Coerced sorted dataset
Out[17]: 100

In [18]: mt2.GT.show()
2018-12-03 11:37:57 Hail: INFO: Coerced sorted dataset
+---------------+------------+-----+------+
| locus | alleles | s | GT |
+---------------+------------+-----+------+
| locus<GRCh37> | array<str> | str | call |
+---------------+------------+-----+------+
| 1:1 | ["A","C"] | "0" | 0/0 |
| 1:1 | ["A","C"] | "1" | 0/0 |
| 1:1 | ["A","C"] | "2" | 0/0 |
| 1:1 | ["A","C"] | "3" | 0/0 |
| 1:1 | ["A","C"] | "4" | 0/0 |
| 1:1 | ["A","C"] | "5" | 0/0 |
| 1:1 | ["A","C"] | "6" | 0/0 |
| 1:1 | ["A","C"] | "7" | 0/0 |
| 1:1 | ["A","C"] | "8" | 0/0 |
| 1:1 | ["A","C"] | "9" | 0/0 |
+---------------+------------+-----+------+
showing top 10 rows

```

Not sure why we coerce the sorted dataset twice for count\_cols.

---

<div class="post-metadata">

### Author: ![tpoterba](https://yyz2.discourse-cdn.com/flex036/user_avatar/discuss.hail.is/tpoterba/32/61_2.png) [@tpoterba](https://discuss.hail.is/u/tpoterba)
#### Post date: [December 3, 2018, 4:51pm UTC](https://discuss.hail.is/t/errors-after-iterative-calls-to-union-cols/790/4 "2018-12-03T16:51:11Z")

</div>

[for dev use]

> Not sure why we coerce the sorted dataset twice for count\_cols.

The first one is key\_rows\_by in BN, the second is the MatrixColsTable in countCols
