I have a dataset of genes I am trying to collect data on from public databases, to use as features in machine learning. I am trying to take some features from UCSC genome browser (e.g. number of CpG islands per gene, number of DNase clusters per gene, regulatory enrichment scores etc.) however I am not sure how to control for bias where a gene that is larger in length - and so will then have more CpG islands or higher regulatory enrichment scores simply due to gene length.

Is there a way to correct for gene length when taking/condensing variant data to individual genes?


For reference, my machine learning model aims to predict whether a gene is the most likely to be causal for a disease (out of all the genes given to the model). The model will score the genes as a regression classification between 0 to 1 (0 being least likely to cause disease and 1 being most likely to cause disease). I plan to later further investigate the genes with the highest scores.

The model uses a variety of multi-omic features (e.g. GTEx gene expression the genes have for many tissues, GWAScatalog data, gene intolerance scores, protein-protein interaction data, drug interaction data, phenotypic scores etc.). However, I am missing epigenetic data to describe my genes so I've been looking to collect based on UCSC's variant data (CpG islands, histone modifications, DNase clusters) - however this leads to my gene length problem when I am trying to reliably take data from the variant level.

I've been plotting my features and gene length, and seen that the UCSC epigenetic data does correlate with having a larger gene length if there is a higher count of regulatory sites (0.8 r2 for some), and so this is what I'm looking to correct.

  • $\begingroup$ Does regulatory enrichment score correspond to the number of regulatory SNPs / other regulatory elements per gene? $\endgroup$
    – PPK
    Aug 4, 2020 at 10:25
  • $\begingroup$ I believe so, a score is given per variant in the public database. I am then finding all those scores per gene in my gene list and usually just picking the highest score to be the value I use for that gene, only problem being a longer gene will have more opportunity for variants and more scores. Although, I am more focused on just creating counts of different regulatory sites per gene (e.g. number of CpG islands) and using that count as a feature input for machine learning (so I can not use the scores at all) but I am still stuck with the gene length issue. $\endgroup$
    – DN1
    Aug 4, 2020 at 12:55
  • $\begingroup$ What do you want to do with your machine learning? (what's correct depends on what are your goals) Have you looked at the way cqn package at Bioconductor corrects for gene length? $\endgroup$
    – llrs
    Aug 4, 2020 at 12:55
  • $\begingroup$ I am looking to use machine learning to predict whether a gene is most likely to cause a disease (model scores the gene with a 0-1 likelihood, predicting that score based on multiomic features I collect to describe the gene). I have looked at cqn but its starting point is RNA-seq data which I don't have. I only have gene symbols and the variants/SNPs within those genes that have come from GWAS to associate to the disease. Thank you for your comments, I'll add this info into my queston for clarity. $\endgroup$
    – DN1
    Aug 4, 2020 at 12:58

1 Answer 1


Its very easy, just let the ML sort this out for you and that is its advantage, You're thinking of GLM style calculation where you pre-screen the data with bivariate plots, where there needs to be nice Q-Q plots and low residual.

For ML simply include the gene length as one of your parameters along with CpG etc ... and the ML regression analysis SVC, lasso, ridge, random forest will figure the relationship out between gene length and CpG. You do zero, the ML does everything, hence from a statistical point of view purists object because you don't know the relationshiop the ML has deduced between the variables, but you will get regression weights for non-DNN stuff, which will give you some idea of the impact of length.

There is the issue of transformations and that can be complicated, but I'd try untransformed data first. The only disadvantage of this approach is the user will have to input the gene size when they want to check out your training algorithm.

  • 1
    $\begingroup$ Thank you for this, this is ideal for me. I'd give you more than 100 points if I could (when it lets me I will award the 100 points). I have a final question for you if possible: when talking about this it was originally suggested I look into statistical tests (I have a biology background so have been unsure) like Fisher’s Exact Test or Chi2. I didn't find anything I thought would easily fit my data - is the GLM style calculation you first mention the only type of statistical testing I could do for this? I can't find any commonly used tests which don't start with RNA-seq data. $\endgroup$
    – DN1
    Aug 4, 2020 at 13:24
  • $\begingroup$ My advice is to stick with ML but do this well, and you don't need a 'fit' anything - you definitely don't need to worry about GLM (ML is more trendy in any case). Geting the ML working to a good standard can be tricky and I'll discuss the issues above shortly. $\endgroup$
    – M__
    Aug 4, 2020 at 16:23
  • $\begingroup$ Thank you for this, I will stick with ML. Currently I use xgboost and view the feature importance, is there a way for me to prove xgboost is regulating how it handles certain features in relation to the gene length feature? I'm trying to finding examples of this but potentially this is a new question for me to ask/research. $\endgroup$
    – DN1
    Aug 6, 2020 at 8:30
  • $\begingroup$ I thought xgboost just handles gradient descent for the ML regression. Are you integrating this with sci-kit learn? Using different methods is important here and further discussion is a good idea. $\endgroup$
    – M__
    Aug 6, 2020 at 16:45
  • $\begingroup$ I use xgboost integrated with scikit-learn, would it be better done with xgboost functions only? I also do an initial benchmarking against xgboost, random forest, GBM, decision tree, k-nearest neighbour and logistic regression (elastic net and lasso). Although xgboost is always the top performer (with 10-fold nested cross validation) when I've tried different things (like different feature selection methods or hyperparameter tuning methods). $\endgroup$
    – DN1
    Aug 7, 2020 at 9:39

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.