tensorQTL interaction issue

Question

[Solal]

Hi,

I hope that the question has its place on this board.

I am trying to run the mapping of cis eQTL using an interaction term with tensorQTL. I am running into the following index error when running the software:

>>> cis.map_nominal(genotype_df, variant_df,
...                 phenotype_df,
...                 phenotype_pos_df,
...                 prefix = prefix,
...                 covariates_df=covariates_df,
...                 interaction_df=interaction_s,
...                 maf_threshold_interaction=0.05,
...                 run_eigenmt=True,
...                 output_dir=output,
...                 write_top=True,
...                 write_stats=True)
cis-QTL mapping: nominal associations for all variant-phenotype pairs
  * 569 samples
  * 24665 phenotypes
  * 68 covariates
  * 69763935 variants
  * including 1 interaction term(s)
    * using 0.05 MAF threshold
  * cis-window: ±1,000,000
  * checking phenotypes: 24665/24665
  * Computing associations
    Mapping chromosome 1
    processing phenotype 137/24665Traceback (most recent call last):
  File "/QRISdata/Q5776/software/tensorqtl/tensorqtl/core.py", line 197, in calculate_interaction_nominal
    Xinv = torch.matmul(torch.transpose(X_t, 1, 2), X_t).inverse() # ng x (1+2*ni) x (1+2*ni)
torch._C._LinAlgError: linalg.inv: (Batch element 2704): The diagonal element 3 is zero, the inversion could not be completed because the input matrix is singular.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/QRISdata/Q5776/software/tensorqtl/tensorqtl/cis.py", line 294, in map_nominal
    res = calculate_interaction_nominal(genotypes_t, phenotype_t.unsqueeze(0), interaction_t,
  File "/QRISdata/Q5776/software/tensorqtl/tensorqtl/core.py", line 200, in calculate_interaction_nominal
    i = int(re.findall('For batch (\d+)', str(e))[0])
IndexError: list index out of range

However, when I run the same command but only on 50 genes, the code runs properly.

>>> cis.map_nominal(genotype_df, variant_df,
...                 phenotype_df.iloc[1:50],
...                 phenotype_pos_df.iloc[1:50],
...                 prefix = prefix,
...                 covariates_df=covariates_df,
...                 interaction_df=interaction_s,
...                 maf_threshold_interaction=0.05,
...                 run_eigenmt=True,
...                 output_dir=output,
...                 write_top=True,
...                 write_stats=True)
cis-QTL mapping: nominal associations for all variant-phenotype pairs
  * 569 samples
  * 49 phenotypes
  * 68 covariates
  * 69763935 variants
  * including 1 interaction term(s)
    * using 0.05 MAF threshold
  * cis-window: ±1,000,000
  * checking phenotypes: 49/49
  * Computing associations
    Mapping chromosome 1
    processing phenotype 49/49
    time elapsed: 0.44 min
    * writing output
done.

The interaction factor is a categorical variable (0/1), I am not sure why the issue is being raised and how to fix

I would appreciate any help/direction to go forward.

Thanks

Answer 1

[Solal]

Writing a solution for anyone interested. The issue stems from a dosage error, some variant contains only heterozygous individuals and no homozygous ones needed for QTL calculation.

I used plink2 to calculate genotypes for all variant using the following command:

plink2  --pfile PATH \
     --geno-counts  \
     --out OUT

I then removed variant that contains less than 10 Homozygous invidivuals (HOM_REF_CT < 10) using plink then was able to run tensorQTL interaction mode.

Hope it helps.