<p>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.</p>
<p>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.</p>
<p>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.</p>

CC BY-SA 4.0
2020-08-04T13:24:23.447Z
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.
DN1
CC BY-SA 4.0
2020-08-04T16:23:31.213Z
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.
M__
CC BY-SA 4.0
2020-08-06T08:30:37.347Z
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.
DN1
CC BY-SA 4.0
2020-08-06T16:45:50.247Z
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.
M__
CC BY-SA 4.0
2020-08-07T09:39:39.727Z
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).
DN1
CC BY-SA 4.0
2020-08-07T13:22:44.107Z
No using sci-kit learn is essential.. What about SVM? What is your accuracy? It is unusual when a single approach easily outperforms everything else in my experience.
M__
CC BY-SA 4.0
2020-08-07T15:54:01.560Z
I've got an r2 of 0.8 with xgboost and 0.75-0.6 r2 with all the other treebased methods. k-nearest is around 0.6 and the logistic regressions are ~0.2 r2. I did a while ago try SVM but found it had extremely low r2 no matter what I did (but I plan to revisit this as I have new/bigger data and hyperparameter tuning). Also I look at other metrics as I know r2 isn't the best, like the mean square error (from 0.16 - 0.18 for the tree-based methods for example). Thanks again for your help as I know this is moving away from the question, regulating gene length has been one of my last puzzle pieces
DN1
CC BY-SA 4.0
2020-08-08T04:49:10.823Z
I'll have another look at your comment later, but I think that r2 is not 'accuracy'. Accuracy is the number of correct positives over the total number of predictions performed on the 'test' split of the data. An accuracy of 0.8 for a biological/genetic system is good enough.and if you have a method to achieve this I think most would suggest this is your result. There may not be universal consensus on this though because some argue that the 'consensus' of all methods is the result that matters (unless it is a DNN).
M__
CC BY-SA 4.0
2020-08-10T09:27:53.533Z
Thank you I'll look into what I can do with getting an accuracy. My only concern is with the regression say the model predicts a gene is scored at 0.801 when the true training score was 0.8 exactly - do I count that as a correct positive or would the model need to be exactly aligned with the real score?
DN1
CC BY-SA 4.0
2020-08-10T12:45:36.810Z
The 0.8 has the look and feel of an accurary score. I wouldn't worry about the 0.801 issue. Accuracy is the primary index used in ML
M__
CC BY-SA 4.0
2020-08-04T13:08:23.503Z
2020-08-04T13:08:23.503Z
M__
<p>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.</p>
<p>Is there a way to correct for gene length when taking/condensing variant data to individual genes?</p>
<p>Edit:</p>
<p>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.</p>
<p>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.</p>
<p>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.</p>

CC BY-SA 4.0
2020-08-04T10:25:18.660Z
Does regulatory enrichment score correspond to the number of regulatory SNPs / other regulatory elements per gene?
PPK
CC BY-SA 4.0
2020-08-04T12:55:41.030Z
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.
DN1
CC BY-SA 4.0
2020-08-04T12:55:51.613Z
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?
llrs
CC BY-SA 4.0
2020-08-04T12:58:42.780Z
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.
DN1
<gene><r><public-databases><genomics><variants>
How to compact variant data to their genes without bias?
# How to compact variant data to their genes without bias?

(Asked by: username_0 on 2020-07-24)

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?

Edit:

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.

**Tags:** gene, r, public-databases, genomics, variants
> **Comments:**
> **username_1**: Does regulatory enrichment score correspond to the number of regulatory SNPs / other regulatory elements per gene?  
> **username_0**: 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.  
> **username_2**: 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?  
> **username_0**: 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.  

---

## 1 answers

### Answer by username_3 on 2020-08-04. Score: 2 | Accepted answer

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.

> **Comments:**
> **username_0**: 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.  
> **username_3**: 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.  
> **username_0**: 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.  
> **username_3**: 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.  
> **username_0**: 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).  
> **username_3**: No using sci-kit learn is essential.. What about SVM? What is your accuracy? It is unusual when a single approach easily outperforms everything else in my experience.  
> **username_0**: I've got an r2 of 0.8 with xgboost and 0.75-0.6 r2 with all the other treebased methods. k-nearest is around 0.6 and the logistic regressions are ~0.2 r2. I did a while ago try SVM but found it had extremely low r2 no matter what I did (but I plan to revisit this as I have new/bigger data and hyperparameter tuning). Also I look at other metrics as I know r2 isn't the best, like the mean square error (from 0.16 - 0.18 for the tree-based methods for example). Thanks again for your help as I know this is moving away from the question, regulating gene length has been one of my last puzzle pieces  
> **username_3**: I'll have another look at your comment later, but I think that r2 is not 'accuracy'. Accuracy is the number of correct positives over the total number of predictions performed on the 'test' split of the data. An accuracy of 0.8 for a biological/genetic system is good enough.and if you have a method to achieve this I think most would suggest this is your result. There may not be universal consensus on this though because some argue that the 'consensus' of all methods is the result that matters (unless it is a DNN).  
> **username_0**: Thank you I'll look into what I can do with getting an accuracy. My only concern is with the regression say the model predicts a gene is scored at 0.801 when the true training score was 0.8 exactly - do I count that as a correct positive or would the model need to be exactly aligned with the real score?  
> **username_3**: The 0.8 has the look and feel of an accurary score. I wouldn't worry about the 0.801 issue. Accuracy is the primary index used in ML  

---

