Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fixing summary calls in simulating.ipynb #239

Merged
merged 1 commit into from
Mar 16, 2016
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
82 changes: 42 additions & 40 deletions documentation_builder/simulating.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,7 @@
{
"data": {
"text/plain": [
"<Solution 0.87 at 0x7f1bb8789550>"
"<Solution 0.87 at 0x10ddd0080>"
]
},
"execution_count": 2,
Expand Down Expand Up @@ -106,7 +106,7 @@
{
"data": {
"text/plain": [
"0.8739215069684909"
"0.8739215069684305"
]
},
"execution_count": 4,
Expand Down Expand Up @@ -172,8 +172,8 @@
" % FLUX RXN ID REACTION \n",
" 41.6% 16 GAPD g3p_c + nad_c + pi_c <=> 13dpg_c + h_c + nadh_c\n",
" 24.1% 9.3 PDH coa_c + nad_c + pyr_c --> accoa_c + co2_c + nadh_c\n",
" 13.1% 5.1 MDH mal__L_c + nad_c <=> h_c + nadh_c + oaa_c\n",
" 13.1% 5.1 AKGDH akg_c + coa_c + nad_c --> co2_c + nadh_c + succoa_c\n",
" 13.1% 5.1 MDH mal__L_c + nad_c <=> h_c + nadh_c + oaa_c\n",
" 8.0% 3.1 Bioma... 1.496 3pg_c + 3.7478 accoa_c + 59.81 atp_c + 0.36...\n",
"\n",
"CONSUMING REACTIONS -- Nicotinamide adenine dinucleotide - reduced\n",
Expand Down Expand Up @@ -263,7 +263,7 @@
{
"data": {
"text/plain": [
"{<Reaction Biomass_Ecoli_core at 0x7f1b92f03470>: 1.0}"
"{<Reaction Biomass_Ecoli_core at 0x116510828>: 1.0}"
]
},
"execution_count": 9,
Expand Down Expand Up @@ -292,7 +292,7 @@
{
"data": {
"text/plain": [
"{<Reaction ATPM at 0x7f1b92f036d8>: 1}"
"{<Reaction ATPM at 0x1165107b8>: 1}"
]
},
"execution_count": 10,
Expand Down Expand Up @@ -320,7 +320,7 @@
{
"data": {
"text/plain": [
"175.00000000002336"
"174.99999999999997"
]
},
"execution_count": 11,
Expand Down Expand Up @@ -349,7 +349,7 @@
{
"data": {
"text/plain": [
"{<Reaction Biomass_Ecoli_core at 0x7f1b92f03470>: 1.0}"
"{<Reaction Biomass_Ecoli_core at 0x116510828>: 1.0}"
]
},
"execution_count": 12,
Expand Down Expand Up @@ -396,17 +396,17 @@
" <tr>\n",
" <th>ACALD</th>\n",
" <td>0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" <tr>\n",
" <th>ACALDt</th>\n",
" <td>0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" <tr>\n",
" <th>ACKr</th>\n",
" <td>-0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" <tr>\n",
" <th>ACONTa</th>\n",
Expand All @@ -421,11 +421,11 @@
" <tr>\n",
" <th>ACt2r</th>\n",
" <td>0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" <tr>\n",
" <th>ADK1</th>\n",
" <td>0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" <tr>\n",
Expand All @@ -436,12 +436,12 @@
" <tr>\n",
" <th>AKGt2r</th>\n",
" <td>0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" <tr>\n",
" <th>ALCD2x</th>\n",
" <td>0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" <tr>\n",
" <th>ATPM</th>\n",
Expand Down Expand Up @@ -475,8 +475,8 @@
" </tr>\n",
" <tr>\n",
" <th>D_LACt2</th>\n",
" <td>-0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" <tr>\n",
" <th>ENO</th>\n",
Expand All @@ -486,11 +486,11 @@
" <tr>\n",
" <th>ETOHt2r</th>\n",
" <td>0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" <tr>\n",
" <th>EX_ac_e</th>\n",
" <td>0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" </tr>\n",
" </tbody>\n",
Expand All @@ -499,26 +499,26 @@
],
"text/plain": [
" maximum minimum\n",
"ACALD 0.00000 -0.00000\n",
"ACALDt 0.00000 -0.00000\n",
"ACKr -0.00000 -0.00000\n",
"ACALD 0.00000 0.00000\n",
"ACALDt -0.00000 0.00000\n",
"ACKr -0.00000 0.00000\n",
"ACONTa 6.00725 6.00725\n",
"ACONTb 6.00725 6.00725\n",
"ACt2r 0.00000 -0.00000\n",
"ADK1 0.00000 0.00000\n",
"ACt2r 0.00000 0.00000\n",
"ADK1 -0.00000 0.00000\n",
"AKGDH 5.06438 5.06438\n",
"AKGt2r 0.00000 -0.00000\n",
"ALCD2x 0.00000 -0.00000\n",
"AKGt2r 0.00000 0.00000\n",
"ALCD2x 0.00000 0.00000\n",
"ATPM 8.39000 8.39000\n",
"ATPS4r 45.51401 45.51401\n",
"Biomass_Ecoli_core 0.87392 0.87392\n",
"CO2t -22.80983 -22.80983\n",
"CS 6.00725 6.00725\n",
"CYTBD 43.59899 43.59899\n",
"D_LACt2 -0.00000 -0.00000\n",
"D_LACt2 0.00000 0.00000\n",
"ENO 14.71614 14.71614\n",
"ETOHt2r 0.00000 -0.00000\n",
"EX_ac_e 0.00000 0.00000"
"ETOHt2r 0.00000 0.00000\n",
"EX_ac_e -0.00000 0.00000"
]
},
"execution_count": 13,
Expand Down Expand Up @@ -561,12 +561,12 @@
" <tbody>\n",
" <tr>\n",
" <th>ACALD</th>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" <td>-2.54237</td>\n",
" </tr>\n",
" <tr>\n",
" <th>ACALDt</th>\n",
" <td>0.00000</td>\n",
" <td>-0.00000</td>\n",
" <td>-2.54237</td>\n",
" </tr>\n",
" <tr>\n",
Expand All @@ -586,7 +586,7 @@
" </tr>\n",
" <tr>\n",
" <th>ACt2r</th>\n",
" <td>-0.00000</td>\n",
" <td>0.00000</td>\n",
" <td>-3.81356</td>\n",
" </tr>\n",
" <tr>\n",
Expand Down Expand Up @@ -665,12 +665,12 @@
],
"text/plain": [
" maximum minimum\n",
"ACALD -0.00000 -2.54237\n",
"ACALDt 0.00000 -2.54237\n",
"ACALD 0.00000 -2.54237\n",
"ACALDt -0.00000 -2.54237\n",
"ACKr -0.00000 -3.81356\n",
"ACONTa 8.89452 0.84859\n",
"ACONTb 8.89452 0.84859\n",
"ACt2r -0.00000 -3.81356\n",
"ACt2r 0.00000 -3.81356\n",
"ADK1 17.16100 0.00000\n",
"AKGDH 8.04593 0.00000\n",
"AKGt2r 0.00000 -1.43008\n",
Expand Down Expand Up @@ -718,7 +718,7 @@
"output_type": "stream",
"text": [
"IN FLUXES OUT FLUXES OBJECTIVES \n",
"o2_e -21.80 ± 1.91 h2o_e 27.86 ± 2.86 Biomass_Ecoli_core 0.000\n",
"o2_e -21.80 ± 1.91 h2o_e 27.86 ± 2.86 Biomass_Ecoli_core 0.874\n",
"glc__D_e -9.76 ± 0.24 co2_e 21.81 ± 2.86 \n",
"nh4_e -4.84 ± 0.32 h_e 19.51 ± 2.86 \n",
"pi_e -3.13 ± 0.08 for_e 2.86 ± 2.86 \n",
Expand All @@ -734,6 +734,7 @@
}
],
"source": [
"model.optimize()\n",
"model.summary(fva=0.95)"
]
},
Expand All @@ -758,13 +759,14 @@
"PRODUCING REACTIONS -- Pyruvate\n",
"-------------------------------\n",
" % FLUX RXN ID REACTION \n",
" 50.0% 6.13 ± 6.13 PYK adp_c + h_c + pep_c --> atp_c + pyr_c\n",
" 50.0% 9.76 ± 0.24 GLCpts glc__D_e + pep_c --> g6p_c + pyr_c\n",
" 85.0% 9.76 ± 0.24 GLCpts glc__D_e + pep_c --> g6p_c + pyr_c\n",
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Those percentages still seem weird. Should fluxes of ~9 and ~6 not end up in perecentages of 9/15 (60%) and 6/15 (40%)? I think the percentages are calculated from the FBA solution and the values from the FVA (but both should be calculated based on the FVA solution).

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Exactly, the percentages are calculated from the FBA solution. This is probably another item to add to #240.

Its non-trivial to decide how percentages should be assigned in the FVA case. Mainly because the total flux through each node is no longer fixed, and the percentages will necessarily change. Even using the median values might not make sense, since the sum of the median inputs may not equal the sum of the median outputs.

Finding errors in percentages would probably also be too hard, although that would be a cool feature.

Perhaps they should just not be printed to avoid confusion.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Good point here with the median input/output problem. Also the median solutions might be confusing. Statistically speaking they are not necessarily correct since the fluxes are continuous probability distributions on a bounded domain and the formula used here assumes complete symmetry which probably does not hold. To get "real" population estimates for the medians one would have to sample from the flux cone with the fixed objective. In papers I have only seen the ranges reported (at least if they do not use sampling) since FVA gives no information about the location of the mean or median of the flux distribution.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

True, I used half-spans +/- range/2 because its a familiar format. Probably better to make summary just have "MAX" and "MIN" columns when called with FVA.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You could also do center of min, max
On Mar 16, 2016 12:24, "Peter St. John" [email protected] wrote:

In documentation_builder/simulating.ipynb
#239 (comment):

@@ -758,13 +759,14 @@
"PRODUCING REACTIONS -- Pyruvate\n",
"-------------------------------\n",
" % FLUX RXN ID REACTION \n",

  •  " 50.0%   6.13 ± 6.13      PYK                  adp_c + h_c + pep_c --> atp_c + pyr_c\n",
    
  •  " 50.0%   9.76 ± 0.24   GLCpts                     glc__D_e + pep_c --> g6p_c + pyr_c\n",
    
  •  " 85.0%   9.76 ± 0.24   GLCpts                     glc__D_e + pep_c --> g6p_c + pyr_c\n",
    

True, I used half-spans +/- range/2 because its a familiar format.
Probably better to make summary just have "MAX" and "MIN" columns when
called with FVA.


You are receiving this because you commented.
Reply to this email directly or view it on GitHub
https://github.com/opencobra/cobrapy/pull/239/files/8a60f16e28ddb40995a91b00b3b74e188f8eb93b#r56398731

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That's actually whats currently done:

        fva_results = pd.DataFrame(
            flux_variability_analysis(met.model, met.reactions,
                                      fraction_of_optimum=fva)).T
        half_span = (fva_results.maximum - fva_results.minimum) / 2
        median = fva_results.minimum + half_span

Fluxes then get printed according to median ± half_span

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@pstjohn yes, completely agree. Sorry for spamming comments here. Maybe we should open a gitter to discuss stuff like this...

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

No worries about all the comments. It's all relevant to the discussion of
this bug.
On Mar 16, 2016 1:15 PM, "Christian Diener" [email protected]
wrote:

In documentation_builder/simulating.ipynb
#239 (comment):

@@ -758,13 +759,14 @@
"PRODUCING REACTIONS -- Pyruvate\n",
"-------------------------------\n",
" % FLUX RXN ID REACTION \n",

  •  " 50.0%   6.13 ± 6.13      PYK                  adp_c + h_c + pep_c --> atp_c + pyr_c\n",
    
  •  " 50.0%   9.76 ± 0.24   GLCpts                     glc__D_e + pep_c --> g6p_c + pyr_c\n",
    
  •  " 85.0%   9.76 ± 0.24   GLCpts                     glc__D_e + pep_c --> g6p_c + pyr_c\n",
    

@pstjohn https://github.com/pstjohn yes, completely agree. Sorry for
spamming comments here. Maybe we should open a gitter to discuss stuff like
this...


You are receiving this because you commented.
Reply to this email directly or view it on GitHub
https://github.com/opencobra/cobrapy/pull/239/files/8a60f16e28ddb40995a91b00b3b74e188f8eb93b#r56406338

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ok :) I would vote for the min - max format since any center value could be interpreted as the mean/median/most likely flux, all of which are false. Also it could give the impression that taking the solution where each flux takes its center value yields the reported biomass objective, which is also not true...

Also,since we are discussing the summaries, wouldn't it be good to allow calculation of the normal summary with pFBA since the individual FBA fluxes are somewhat arbitrary (pFBA as well but much less).

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So right now summary methods don't trigger an FBA solution, pFBA or otherwise. You can certainly get pFBA solutions in the summary methods, as long as the model is optimized prior to calling the summaries. This is partly by design, as they are supposed to be quick ways of examining a model's current state.

This is muddled somewhat by the addition of the FVA options, which do trigger optimizations.

" 15.0% 6.13 ± 6.13 PYK adp_c + h_c + pep_c --> atp_c + pyr_c\n",
"\n",
"CONSUMING REACTIONS -- Pyruvate\n",
"-------------------------------\n",
" % FLUX RXN ID REACTION \n",
"100.0% 11.34 ± 7.43 PDH coa_c + nad_c + pyr_c --> accoa_c + co2_c + nadh_c\n"
" 78.9% 11.34 ± 7.43 PDH coa_c + nad_c + pyr_c --> accoa_c + co2_c + nadh_c\n",
" 21.1% 0.85 ± 0.02 Bioma... 1.496 3pg_c + 3.7478 accoa_c + 59.81 atp_c + 0.36...\n"
]
}
],
Expand Down Expand Up @@ -817,7 +819,7 @@
{
"data": {
"text/plain": [
"6.072919944699606e-14"
"1.1102230246251565e-16"
]
},
"execution_count": 18,
Expand Down Expand Up @@ -846,7 +848,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.4.3"
"version": "3.4.4"
}
},
"nbformat": 4,
Expand Down