Skip to content

Commit 6e60131

Browse files
author
Scott Campit
committed
First commit
0 parents  commit 6e60131

25 files changed

+1213
-0
lines changed

CODE_OF_CONDUCT.md

Lines changed: 34 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,34 @@
1+
Contributor Covenant Code of Conduct
2+
====================================
3+
4+
## Our Pledge
5+
In the interest of fostering an open and welcoming environment, we as contributors and maintainers pledge to making participation in our project and our community a harassment-free experience for everyone, regardless of age, body size, disability, ethnicity, gender identity and expression, level of experience, nationality, personal appearance, race, religion, or sexual identity and orientation.
6+
7+
## Our Standards
8+
Examples of behavior that contributes to creating a positive environment include:
9+
* Using welcoming and inclusive language
10+
* Being respectful of differing viewpoints and experiences
11+
* Gracefully accepting constructive criticism
12+
* Focusing on what is best for the community
13+
* Showing empathy towards other community members
14+
15+
Examples of unacceptable behavior by participants include:
16+
* The use of sexualized language or imagery and unwelcome sexual attention or advances
17+
* Trolling, insulting/derogatory comments, and personal or political attacks
18+
* Public or private harassment
19+
* Publishing others' private information, such as a physical or electronic address, without explicit permission
20+
* Other conduct which could reasonably be considered inappropriate in a professional setting
21+
22+
## Our Responsibilities
23+
Project contributors are responsible for clarifying the standards of acceptable behavior and are expected to take appropriate and fair corrective action in response to any instances of unacceptable behavior.
24+
Project contributors have the right and responsibility to remove, edit, or reject comments, commits, code, wiki edits, issues, and other contributions that are not aligned to this Code of Conduct, or to ban temporarily or permanently any contributor for other behaviors that they deem inappropriate, threatening, offensive, or harmful.
25+
26+
## Scope
27+
This Code of Conduct applies both within project spaces and in public spaces when an individual is representing the project or its community. Examples of representing a project or community include using an official project e-mail address, posting via an official social media account, or acting as an appointed representative at an online or offline event. Representation of a project may be further defined and clarified by project contributors.
28+
29+
## Enforcement
30+
If a contributor engages in harassing behaviour, the project organizers may take any action they deem appropriate, including warning the offender or expelling them from online forums, online project resources, face-to-face meetings, or any other project-related activity or resource.
31+
If you are being harassed, notice that someone else is being harassed, or have any other concerns, please contact a member of the project team immediately. Instances of abusive, harassing, or otherwise unacceptable behavior may be reported by contacting the project team. All complaints will be reviewed and investigated and will result in a response that is deemed necessary and appropriate to the circumstances. The project team is obligated to maintain confidentiality with regard to the reporter of an incident. Further details of specific enforcement policies may be posted separately.
32+
33+
## Attribution
34+
Portions of this Code of Conduct were adapted from Electron's [Contributor Covenant Code of Conduct](https://github.com/electron/electron/blob/master/CODE_OF_CONDUCT.md), which itself was adapted from the [Contributor Covenant](http://contributor-covenant.org/version/1/4), version 1.4.

CONTRIBUTING.md

Lines changed: 11 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,11 @@
1+
# Guidelines for contributing to this project
2+
Any constructive contributions – bug reports, pull requests (code or documentation), suggestions for improvements, and more – are welcome.
3+
4+
## Conduct
5+
Everyone is asked to read and respect the [code of conduct](CODE_OF_CONDUCT.md) before participating in this project.
6+
7+
## Coordinating work
8+
A quick way to find out what is currently in the near-term plans for this project is to look at the [GitHub issue tracker](https://github.com/mhucka/readmine/issues), but the possibilities are not limited to what you see there – if you have ideas for new features and enhancements, please feel free to write them up as a new issue or contact the developers directly!
9+
10+
## Submitting contributions
11+
Please feel free to contact the author directly, or even better, jump right in and use the standard GitHub approach of forking the repo and creating a pull request. When committing code changes and submitting pull requests, please write a clear log message for your commits.

README.md

Lines changed: 52 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,52 @@
1+
# Inferring metabolic fluxes from time-course metabolomics data using Dynamic Flux Analysis (DFA)
2+
Metabolism is highly dynamic, and is tied to essential cellular processes such as cellular differentiation and proliferation. However, it is challenging to measure the activity of every reaction in a cell. Time-course metabolomics is a useful measure for tracing metabolite dynamics, but most approaches that utilize this data set examines a small part of the metabolic network.
3+
4+
**Dynamic flux analysis (DFA)** addresses this issue by fitting time-course metabolomics with genome-scale metabolic models to infer metabolic activity throughout the entire metabolic network.
5+
6+
*Advantages of DFA include*:
7+
* **Flexibility of data type:** DFA can use both absolute and relative metabolomics values.
8+
* **Handling large and small data:** DFA can use as little as two time points, while other models require several time points to fit statistical or genome-scale model.
9+
* **Ease of parameter tuning:** DFA models are easily tuned by a single parameter **kappa**, which is a weight proportional to the optimization strength of the objective function.
10+
11+
## Installation
12+
DFA was implemented in [MATLAB](https://www.mathworks.com/products/matlab.html) (recommended version: 2018+), and requires the [Gurobi Mathematical Programming Solver](https://www.gurobi.com/) (recommended version: 2018+). The file `dfa.m` or the livescript `dfa.mlx` in this repository are used to run DFA.
13+
14+
## Usage
15+
The function `dfa` has the following syntax:
16+
```MATLAB
17+
[dfaModel, solution] = DFA(model, metabolomics, params)
18+
```
19+
20+
**INPUTS:**
21+
* `model`: A `structure` of the genome-scale metabolic model in COBRA format.
22+
* `metabolomics`: A `structure` of the metabolomics data. It requires the following `fields`:
23+
* `data`: A numerical array of the time course metabolomics data. Rows correspond to metabolites, columns correspond to individual time points.
24+
* `position`: A numerical array mapping the metabolite to compartments in the metabolic model. Rows correspond to metabolites, columns correspond to individual compartments.
25+
* `params (optional)`: A `structure` containing parameters and hyperparameters. It takes the following `fields`:
26+
* `kappa`: A scalar denoting the weight to penalize the metabolomics data. The default value is 1.
27+
* `kappa2`: A scalar denoting the weights to minimize the sum of fluxes using parsimonious Flux Balance Analysis (pFBA). The default value is 1E-6.
28+
* `norm`: A string denoting which function to use for normalizing the flux activity coefficients. There are three possible arguments. The default argument is `None`:
29+
* `None`: No normalization for the flux activity coefficients
30+
* `MAV`: Normalizing the flux activity coefficients from the maximum absolute value
31+
* `Quantile`: Quantile normalization of the flux activity coefficients
32+
* `pfba`: A Boolean that determines whether to run pFBA on or off.
33+
34+
**OUTPUTS:**
35+
* `dfaModel`: A `structure` of the genome-scale metabolic model in COBRA format fitted with the flux activity coefficients computed from the time-course metabolomics.
36+
* `solution`: A `structure` containing solutions from the Gurobi solver, including the metabolic fluxes.
37+
38+
## Contributing
39+
Contributions are welcome! Please read the contributions guide to get started. Also feel free to submit bugs, feature requests, and pull requests.
40+
41+
Additionally, you can support development for DFA by citing the original manuscript:
42+
43+
Chandrasekaran, Sriram, Jin Zhang, Zhen Sun, Li Zhang, Christian A. Ross, Yu-Chung Huang, John M. Asara, Hu Li, George Q. Daley, and James J. Collins. ["Comprehensive mapping of pluripotent stem cell metabolism using dynamic genome-scale network modeling."](https://www.cell.com/cell-reports/fulltext/S2211-1247(17)31027-6) Cell reports 21, no. 10 (2017): 2965-2977.
44+
45+
## Publications using Dynamic Flux Analysis
46+
**Stem cell metabolism**
47+
1. Chandrasekaran, Sriram, Jin Zhang, Zhen Sun, Li Zhang, Christian A. Ross, Yu-Chung Huang, John M. Asara, Hu Li, George Q. Daley, and James J. Collins. ["Comprehensive mapping of pluripotent stem cell metabolism using dynamic genome-scale network modeling."](https://www.cell.com/cell-reports/fulltext/S2211-1247(17)31027-6) Cell reports 21, no. 10 (2017): 2965-2977.
48+
2. Shen, Fangzhou, Camden Cheek, and Sriram Chandrasekaran. ["Dynamic network modeling of stem cell metabolism."](https://link.springer.com/protocol/10.1007%2F978-1-4939-9224-9_14) In Computational stem cell biology, pp. 305-320. Humana, New York, NY, 2019.
49+
50+
**Cancer metabolism**
51+
1. Nelson, Barbara S., Lin Lin, Daniel M. Kremer, Cristovão M. Sousa, Cecilia Cotta-Ramusino, Amy Myers, Johanna Ramos et al. ["Tissue of origin dictates GOT1 dependence and confers synthetic lethality to radiotherapy."](https://cancerandmetabolism.biomedcentral.com/articles/10.1186/s40170-019-0202-2) Cancer & Metabolism 8, no. 1 (2020): 1-16.
52+
2. Campit, Scott, and Sriram Chandrasekaran. ["Inferring Metabolic Flux from Time-Course Metabolomics."](https://link.springer.com/protocol/10.1007%2F978-1-0716-0159-4_13) In Metabolic Flux Analysis in Eukaryotic Cells, pp. 299-313. Humana, New York, NY, 2020.

data/tutorial.xlsx

210 KB
Binary file not shown.

model/model_human_duarte.mat

328 KB
Binary file not shown.

src/DFA.m

Lines changed: 155 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,155 @@
1+
function [dynamicModel, solution] = DFA(model, metabolomics, params)
2+
%% DFA Dynamic flux analysis for time-course metabolomics data
3+
% |dfa| computes metabolic fluxes using time-course metabolomics data.
4+
%
5+
% *INPUTS*
6+
%%
7+
% # |model|: A geneome-scale metabolic model in COBRA format
8+
% # |file|: An excel spreadsheet containing the metabolomics
9+
% dataset.
10+
% # |sheetname|: The sheet name in the Excel spreadsheet. The default
11+
% sheet is 'Sheet1'
12+
% # |kappa|: The relative weight for the consistency of the flux
13+
% distribution compared to the weight for maximizing the objective function. The
14+
% default value is 1
15+
% # |kappa2|: The relative weight for minimizing the sum of absolute
16+
% flux values of other reactions. The default value is 1E-3
17+
%%
18+
% *OUTPUTS*
19+
%%
20+
% # |dynamicModel|: A metabolic model with pseudoreactions that correspond
21+
% to the time-course metabolomics data
22+
% # |solution|: The growth rate objective value from single gene KO
23+
24+
% Default parameters
25+
if (~exist('params.kappa', 'var')) || (isempty(params.kappa))
26+
params.kappa = 1;
27+
end
28+
if (~exist('params.kappa2', 'var')) || (isempty(params.kappa2))
29+
params.kappa2 = 1E-3;
30+
end
31+
if (~exist('params.norm', 'var')) || (isempty(params.norm))
32+
params.norm = 'None';
33+
end
34+
35+
% Initialize Gurobi model parameters
36+
tmpModel = model;
37+
tmpModel.A = model.S;
38+
tmpModel.obj = model.c;
39+
tmpModel.rhs = model.b;
40+
tmpModel.sense = repmat('=', [size(model.S, 1), 1]);
41+
tmpModel.lb = model.lb;
42+
tmpModel.ub = model.ub;
43+
tmpModel.vtype = repmat('C', [size(model.S, 2), 1]);
44+
tmpModel.modelsense = 'max';
45+
tmpModel.original_size = length(tmpModel.rxns);
46+
47+
% Compute the linear flux activity coefficient
48+
slope = zeros(size(metabolomics.data, 1), 1);
49+
for m = 1:size(metabolomics.data, 1)
50+
weight = polyfit(1:size(metabolomics.data, 2), metabolomics.data(m, :), 1);
51+
slope(m) = weight(1) / weight(2);
52+
end
53+
54+
% Normalize the flux activity coefficients
55+
switch params.norm
56+
case 'None'
57+
slope = slope;
58+
case 'MAV'
59+
slope(:, 1) = slope(:, 1) ./ max(abs(slope(:, 1)));
60+
case 'Quantile'
61+
slope(:, 1) = slope(:, 1) ./ quantilenorm(abs(slope(:, 1)), 'MEDIAN', true);
62+
slope(slope(:, 1) > 1, 1) = log10(slope(slope(:, 1) > 1, 1)) + 1;
63+
slope(slope(:, 1) < -1, 1) = -log10(-1 * slope(slope(:, 1) < -1, 1)) - 1;
64+
end
65+
66+
% Create the right hand side and constrain the model using the flux
67+
% activity coefficients
68+
constrainedModel = tmpModel;
69+
epsilon = zeros(size(slope(:, 1)));
70+
for m = 1:size(metabolomics.positions, 1)
71+
metabolite = (metabolomics.positions ~= 0);
72+
if ~isempty(metabolite)
73+
positions = metabolomics.positions(m, :);
74+
positions(positions == 0) = [];
75+
epsilon(m) = slope(m);
76+
77+
alpha = size(constrainedModel.A, 2) + 1;
78+
beta = alpha + 1;
79+
gamma = size(constrainedModel.A, 1) + 1;
80+
constrainedModel.rhs(positions) = epsilon(m);
81+
constrainedModel.rhs(isnan(constrainedModel.rhs)) = 0;
82+
83+
% metabolite constraint = tolerance - alpha + beta
84+
constrainedModel.A(positions, alpha) = 1;
85+
constrainedModel.A(positions, beta) = -1;
86+
constrainedModel.vtype(alpha) = 'C';
87+
constrainedModel.lb(beta) = 0;
88+
constrainedModel.lb(alpha) = 0;
89+
constrainedModel.ub(beta) = 1000;
90+
constrainedModel.ub(alpha) = 1000;
91+
constrainedModel.obj(beta) = -1 * params.kappa;
92+
constrainedModel.obj(alpha) = -1 * params.kappa;
93+
constrainedModel.c(alpha) = -1 * params.kappa;
94+
constrainedModel.c(beta) = -1 * params.kappa;
95+
constrainedModel.rxns{alpha} = 'alpha';
96+
constrainedModel.rxns{beta} = 'beta';
97+
constrainedModel.mets{gamma} = 'pmet';
98+
constrainedModel.sense(alpha) = '=';
99+
constrainedModel.sense(beta) = '=';
100+
101+
end
102+
end
103+
104+
% Minimize the sum of absolute flux values for all other reactions
105+
% using kappa2 (pFBA)
106+
for r = 1:length(constrainedModel.rxns)
107+
if constrainedModel.c(r) == 0
108+
109+
% xi + si >= -eps2
110+
% si >= 0
111+
% rho(ri + si)
112+
pmet = size(constrainedModel.A, 1) + 1;
113+
prxn = size(constrainedModel.A, 2) + 1;
114+
constrainedModel.A(pmet, r) = 1;
115+
constrainedModel.A(pmet, prxn) = 1;
116+
constrainedModel.rhs(pmet) = -0;
117+
constrainedModel.sense(pmet) = '>';
118+
constrainedModel.lb(prxn) = 0;
119+
constrainedModel.ub(prxn) = 1000;
120+
constrainedModel.obj(prxn) = -1 * params.kappa2;
121+
constrainedModel.c(prxn) = -1 * params.kappa2;
122+
constrainedModel.rxns{prxn} = 'pseudorxn';
123+
% xi - ri <= eps2
124+
% ri >= 0
125+
pmet = size(constrainedModel.A, 1) + 1;
126+
prxn = size(constrainedModel.A, 2) + 1;
127+
constrainedModel.A(pmet, r) = 1;
128+
constrainedModel.A(pmet, prxn) = -1;
129+
constrainedModel.rhs(pmet) = 0;
130+
constrainedModel.sense(pmet) = '<';
131+
constrainedModel.lb(prxn) = 0;
132+
constrainedModel.ub(prxn) = 1000;
133+
constrainedModel.obj(prxn) = -1 * params.kappa2;
134+
constrainedModel.rxns{prxn} = 'pseudorxn';
135+
end
136+
end
137+
138+
constrainedModel.vtype = repmat('C', size(constrainedModel.A, 2), 1);
139+
dynamicModel = constrainedModel;
140+
141+
% Set params for Gurobi
142+
params2.outputflag = 0;
143+
params2.Threads = 2;
144+
params2.Seed = 314;
145+
params2.NumericFocus = 3;
146+
solution = gurobi(dynamicModel, params2);
147+
148+
try
149+
solution.v = solution.x(1:length(model.rxns));
150+
solution.solverObj = solution.objval;
151+
solution.rhs = slope;
152+
catch ME
153+
warning('Infeasible solution obtained for a metabolic model')
154+
end
155+
end

src/DFA.mlx

5.71 KB
Binary file not shown.

src/EDA/UnitTestForLoop.mlx

3.32 KB
Binary file not shown.

src/EDA/all_CV_Filter.mlx

4 KB
Binary file not shown.

src/EDA/computeReplicateMetrics.mlx

3.98 KB
Binary file not shown.

0 commit comments

Comments
 (0)