[CHNOSZ-commits] r896 - in pkg/CHNOSZ: . R demo inst inst/tinytest man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun May 11 13:20:20 CEST 2025
Author: jedick
Date: 2025-05-11 13:20:20 +0200 (Sun, 11 May 2025)
New Revision: 896
Removed:
pkg/CHNOSZ/demo/protein.equil.R
pkg/CHNOSZ/demo/protein_A.R
pkg/CHNOSZ/demo/protein_Cp.R
Modified:
pkg/CHNOSZ/DESCRIPTION
pkg/CHNOSZ/NAMESPACE
pkg/CHNOSZ/R/examples.R
pkg/CHNOSZ/R/protein.info.R
pkg/CHNOSZ/demo/00Index
pkg/CHNOSZ/inst/NEWS.Rd
pkg/CHNOSZ/inst/TODO
pkg/CHNOSZ/inst/tinytest/test-protein.info.R
pkg/CHNOSZ/man/examples.Rd
pkg/CHNOSZ/man/ionize.aa.Rd
pkg/CHNOSZ/man/protein.info.Rd
Log:
Move protein.equil() to JMDplots
Modified: pkg/CHNOSZ/DESCRIPTION
===================================================================
--- pkg/CHNOSZ/DESCRIPTION 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/DESCRIPTION 2025-05-11 11:20:20 UTC (rev 896)
@@ -1,6 +1,6 @@
-Date: 2025-05-10
+Date: 2025-05-11
Package: CHNOSZ
-Version: 2.1.0-67
+Version: 2.1.0-68
Title: Thermodynamic Calculations and Diagrams for Geochemistry
Authors at R: c(
person("Jeffrey", "Dick", , "j3ffdick at gmail.com", role = c("aut", "cre"),
Modified: pkg/CHNOSZ/NAMESPACE
===================================================================
--- pkg/CHNOSZ/NAMESPACE 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/NAMESPACE 2025-05-11 11:20:20 UTC (rev 896)
@@ -28,7 +28,7 @@
"ratlab",
"EOSregress", "EOScoeffs", "EOSplot", "EOSvar",
# demos
- "protein.equil", "palply",
+ "palply",
"label.plot",
"basis.logact",
"label.figure", "syslab",
Modified: pkg/CHNOSZ/R/examples.R
===================================================================
--- pkg/CHNOSZ/R/examples.R 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/R/examples.R 2025-05-11 11:20:20 UTC (rev 896)
@@ -30,12 +30,12 @@
cat("Time elapsed: ", proc.time() - .ptime, "\n")
}
-demos <- function(which = c("references", "protein.equil", "affinity", "NaCl", "density",
+demos <- function(which = c("references", "affinity", "NaCl", "density",
"ORP", "ionize", "buffer", "protbuff", "glycinate",
"mosaic", "copper", "arsenic", "solubility", "gold", "contour", "sphalerite", "minsol",
"Shh", "saturation", "adenine", "DEW", "lambda", "potassium", "TCA", "aluminum",
"AD", "comproportionation", "Pourbaix", "E_coli", "yttrium", "rank.affinity", "uranyl",
- "sum_S", "MgATP", "protein_Cp", "protein_A", "rubisco_Zc"),
+ "sum_S", "MgATP", "rubisco_Zc"),
save.png = FALSE) {
# Run one or more demos from CHNOSZ with ask = FALSE, and return the value of the last one
out <- NULL
Modified: pkg/CHNOSZ/R/protein.info.R
===================================================================
--- pkg/CHNOSZ/R/protein.info.R 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/R/protein.info.R 2025-05-11 11:20:20 UTC (rev 896)
@@ -6,7 +6,6 @@
# protein.formula: chemical makeup of the indicated proteins
# protein.OBIGT: perform group additivity calculations
# protein.basis: coefficients of basis species in formation reactions of [ionized] proteins [residues]
-# protein.equil: step-by-step example of protein equilibrium calculation
pinfo <- function(protein, organism = NULL, residue = FALSE, regexp = FALSE) {
# Return the `protein` (possibly per residue) for:
@@ -151,122 +150,3 @@
return(sb)
}
-protein.equil <- function(protein, T = 25, loga.protein = 0, digits = 4) {
- out <- character()
- mymessage <- function(...) {
- message(...)
- text <- paste(list(...), collapse = " ")
- out <<- c(out, text)
- }
- # Show the individual steps in calculating metastable equilibrium among proteins
- mymessage("protein.equil: temperature from argument is ", T, " degrees C")
- # Display units
- E_units <- E.units()
- mymessage("protein.equil: energy units is ", E_units)
- TK <- convert(T, "K")
- # Get the amino acid compositions of the proteins
- aa <- pinfo(pinfo(protein))
- # Get some general information about the proteins
- pname <- paste(aa$protein, aa$organism, sep = "_")
- plength <- protein.length(aa)
- # Use thermo()$basis to decide whether to ionize the proteins
- thermo <- get("thermo", CHNOSZ)
- ionize.it <- FALSE
- iword <- "nonionized"
- bmat <- basis.elements()
- if("H+" %in% rownames(bmat)) {
- ionize.it <- TRUE
- iword <- "ionized"
- pH <- -thermo$basis$logact[match("H+", rownames(bmat))]
- mymessage("protein.equil: pH from thermo$basis is ", pH)
- }
- # Tell the user whose [Met] is in thermo$OBIGT
- info.Met <- info(info('[Met]', "aq"))
- mymessage("protein.equil: [Met] is from reference ", info.Met$ref1)
- ## First set of output: show results of calculations for a single protein
- mymessage("protein.equil [1]: first protein is ", pname[1], " with length ", plength[1])
- # Standard Gibbs energies of basis species
- G0basis <- unlist(suppressMessages(subcrt(thermo$basis$ispecies, T = T, property = "G")$out))
- # Coefficients of basis species in formation reactions of proteins
- protbasis <- suppressMessages(protein.basis(aa, T = T))
- # Sum of standard Gibbs energies of basis species in each reaction
- G0basissum <- colSums(t(protbasis) * G0basis)
- # Standard Gibbs energies of nonionized proteins
- G0prot <- unlist(suppressMessages(subcrt(pname, T = T, property = "G")$out))
- # Standard Gibbs energy of formation reaction of nonionized protein, E_units/mol
- G0protform <- G0prot - G0basissum
- mymessage("protein.equil [1]: reaction to form nonionized protein from basis species has G0(", E_units, "/mol) of ", signif(G0protform[1], digits))
- if(ionize.it) {
- # Standard Gibbs energy of ionization of protein, J/mol
- G0ionization <- suppressMessages(ionize.aa(aa, property = "G", T = T, pH = pH))[1, ]
- # Standard Gibbs energy of ionization of protein, E_units/mol
- if(E_units == "cal") G0ionization <- convert(G0ionization, "cal")
- mymessage("protein.equil [1]: ionization reaction of protein has G0(", E_units, "/mol) of ", signif(G0ionization[1], digits))
- # Standard Gibbs energy of formation reaction of ionized protein, E_units/mol
- G0protform <- G0protform + G0ionization
- }
- # Standard Gibbs energy of formation reaction of non/ionized residue equivalents, dimensionless
- # Gas constant
- #if(E_units == "cal") R <- 1.9872 # gas constant, cal K^-1 mol^-1
- #if(E_units == "J") R <- 8.314445 # = 1.9872 * 4.184 J K^-1 mol^-1 20220325
- if(E_units == "J") R <- 8.314463 # https://physics.nist.gov/cgi-bin/cuu/Value?r 20230630
- if(E_units == "cal") R <- 8.314463 / 4.184
- G0res.RT <- G0protform/R/TK/plength
- mymessage("protein.equil [1]: per residue, reaction to form ", iword, " protein from basis species has G0/RT of ", signif(G0res.RT[1], digits))
- # Coefficients of basis species in formation reactions of residues
- resbasis <- suppressMessages(protein.basis(aa, T = T, normalize = TRUE))
- # logQstar and Astar/RT
- logQstar <- colSums(t(resbasis) * - thermo$basis$logact)
- mymessage("protein.equil [1]: per residue, logQstar is ", signif(logQstar[1], digits))
- Astar.RT <- -G0res.RT - log(10)*logQstar
- mymessage("protein.equil [1]: per residue, Astar/RT = -G0/RT - 2.303logQstar is ", signif(Astar.RT[1], digits))
- if(!is.numeric(protein)) mymessage("protein.equil [1]: not comparing calculations with affinity() because 'protein' is not numeric")
- else {
- # For **Astar** we have to set the activities of the proteins to zero, not loga.protein!
- a <- suppressMessages(affinity(iprotein = protein, T = T, loga.protein = 0))
- aAstar.RT <- log(10) * as.numeric(a$values) / plength
- mymessage("check it! per residue, Astar/RT calculated using affinity() is ", signif(aAstar.RT[1], digits))
- if(!isTRUE(all.equal(Astar.RT, aAstar.RT, check.attributes = FALSE)))
- stop("Bug alert! The same value for Astar/RT cannot be calculated manually as by using affinity()")
- }
- if(length(pname) == 1) mymessage("protein.equil [all]: all done... give me more than one protein for equilibrium calculations")
- else {
- ## Next set of output: equilibrium calculations
- mymessage("protein.equil [all]: lengths of all proteins are ", paste(plength, collapse = " "))
- mymessage("protein.equil [all]: Astar/RT of all residue equivalents are ", paste(signif(Astar.RT, digits), collapse = " "))
- expAstar.RT <- exp(Astar.RT)
- sumexpAstar.RT <- sum(expAstar.RT)
- mymessage("protein.equil [all]: sum of exp(Astar/RT) of all residue equivalents is ", signif(sumexpAstar.RT, digits))
- # Boltzmann distribution
- alpha <- expAstar.RT / sumexpAstar.RT
- mymessage("protein.equil [all]: equilibrium degrees of formation (alphas) of residue equivalents are ", paste(signif(alpha, digits), collapse = " "))
- # Check with equilibrate()
- if(is.numeric(protein)) {
- loga.equil.protein <- unlist(suppressMessages(equilibrate(a, normalize = TRUE))$loga.equil)
- # Here we do have to convert from logarithms of activities of proteins to degrees of formation of residue equivalents
- a.equil.residue <- plength*10^loga.equil.protein
- ealpha <- a.equil.residue/sum(a.equil.residue)
- mymessage("check it! alphas of residue equivalents from equilibrate() are ", paste(signif(ealpha, digits), collapse = " "))
- if(!isTRUE(all.equal(alpha, ealpha, check.attributes = FALSE)))
- stop("Bug alert! The same value for alpha cannot be calculated manually as by using equilibrate()")
- }
- # Total activity of residues
- loga.residue <- log10(sum(plength * 10^loga.protein))
- mymessage("protein.equil [all]: for activity of proteins equal to 10^", signif(loga.protein, digits), ", total activity of residues is 10^", signif(loga.residue, digits))
- # Equilibrium activities of residues
- loga.residue.equil <- log10(alpha*10^loga.residue)
- mymessage("protein.equil [all]: log10 equilibrium activities of residue equivalents are ", paste(signif(loga.residue.equil, digits), collapse = " "))
- # Equilibrium activities of proteins
- loga.protein.equil <- log10(10^loga.residue.equil/plength)
- mymessage("protein.equil [all]: log10 equilibrium activities of proteins are ", paste(signif(loga.protein.equil, digits), collapse = " "))
- # Check with equilibrate()
- if(is.numeric(protein)) {
- eloga.protein.equil <- unlist(suppressMessages(equilibrate(a, loga.balance = loga.residue, normalize = TRUE))$loga.equil)
- mymessage("check it! log10 eq'm activities of proteins from equilibrate() are ", paste(signif(eloga.protein.equil, digits), collapse = " "))
- if(!isTRUE(all.equal(loga.protein.equil, eloga.protein.equil, check.attributes = FALSE)))
- stop("Bug alert! The same value for log10 equilibrium activities of proteins cannot be calculated manually as by using equilibrate()")
- }
- }
- return(out)
-}
-
Modified: pkg/CHNOSZ/demo/00Index
===================================================================
--- pkg/CHNOSZ/demo/00Index 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/demo/00Index 2025-05-11 11:20:20 UTC (rev 896)
@@ -1,5 +1,4 @@
references Cross-check the references in refs.csv with the thermodynamic database
-protein.equil Chemical activities of two proteins in metastable equilibrium
affinity Affinities of metabolic reactions and amino acid synthesis
NaCl Equilibrium constant for aqueous NaCl dissociation
density Density of H2O, inverted from IAPWS-95 equations
@@ -35,6 +34,4 @@
yttrium logK.to.OBIGT fits at 800 and 1000 bar and Y speciation in NaCl solution at varying pH
uranyl Total (carbonate|sulfate)-pH diagrams for uranyl species
MgATP Speciation of ATP with H+ and Mg+2
-protein_Cp Compare calculated and experimental Cp of proteins
-protein_A Affinity of ionization of proteins
rubisco_Zc Zc of Rubisco vs optimal growth temperature
Deleted: pkg/CHNOSZ/demo/protein.equil.R
===================================================================
--- pkg/CHNOSZ/demo/protein.equil.R 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/demo/protein.equil.R 2025-05-11 11:20:20 UTC (rev 896)
@@ -1,32 +0,0 @@
-# CHNOSZ/demo/protein.equil.R
-## Steps in calculation of chemical activities of two proteins
-## in metastable equilibrium, after Dick and Shock, 2011
-library(CHNOSZ)
-
-protein <- pinfo(c("CSG_METVO", "CSG_METJA"))
-# Use superseded properties of [Met], [Gly], and [UPBB] (Dick et al., 2006)
-mod.OBIGT("[Met]", G = -35245, H = -59310, S = 40.38)
-mod.OBIGT("[Gly]", G = -6075, H = -5570, S = 17.31)
-mod.OBIGT("[UPBB]", G = -21436, H = -45220, S = 1.62)
-# Set up the basis species to those used in DS11
-basis("CHNOS+")
-# Note this yields logaH2 = -4.657486
-swap.basis("O2", "H2")
-# Demonstrate the steps of the equilibrium calculation
-protein.equil(protein, loga.protein = -3)
-## We can also look at the affinities
-# (Reaction 7, Dick and Shock, 2011)
-# A/2.303RT for protein at unit activity (A-star for the protein)
-a <- affinity(iprotein = protein[1], loga.protein = 0)
-Astar.protein <- a$values[[1]]
-# Divide affinity by protein length (A-star for the residue)
-pl <- protein.length(protein[1])
-Astar.residue <- a$values[[1]]/pl # 0.1893, Eq. 11
-# A/2.303RT per residue corresponding to protein activity of 10^-3
-loga.residue <- log10(pl*10^-3)
-Aref.residue <- Astar.residue - loga.residue # 0.446, after Eq. 16
-# A-star of the residue in natural log units (A/RT)
-log(10) * Astar.residue # 0.4359, after Eq. 23
-
-# Forget about the superseded group properties for whatever comes next
-reset()
Deleted: pkg/CHNOSZ/demo/protein_A.R
===================================================================
--- pkg/CHNOSZ/demo/protein_A.R 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/demo/protein_A.R 2025-05-11 11:20:20 UTC (rev 896)
@@ -1,23 +0,0 @@
-# CHNOSZ/demo/protein_A.R
-# Calculate affinity of ionization of proteins
-# 20250510 demo extracted from anintro.Rmd
-
-library(CHNOSZ)
-
-ip <- pinfo(c("CYC_BOVIN", "LYSC_CHICK", "MYG_PHYCA", "RNAS1_BOVIN"))
-basis("CHNOS+")
-a_ion <- affinity(pH = c(0, 14), iprotein = ip)
-basis("CHNOS")
-a_nonion <- affinity(iprotein = ip)
-plot(c(0, 14), c(50, 300), xlab = "pH", ylab = quote(italic(A/2.303*RT)), type = "n")
-for(i in 1:4) {
- A_ion <- as.numeric(a_ion$values[[i]])
- A_nonion <- as.numeric(a_nonion$values[[i]])
- lines(a_ion$vals[[1]], A_ion - A_nonion, col=i)
-}
-legend("topright", legend = a_ion$species$name,
- col = 1:4, lty = 1, bty = "n", cex = 0.9)
-title("Affinity of ionization of proteins", font.main = 1)
-
-# The affinity is always positive, showing the strong energetic drive for ionization of proteins in aqueous solution.
-# The degrees of ionization of amino and carboxyl groups increase at low and high pH, respectively, giving rise to the U-shaped lines.
Deleted: pkg/CHNOSZ/demo/protein_Cp.R
===================================================================
--- pkg/CHNOSZ/demo/protein_Cp.R 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/demo/protein_Cp.R 2025-05-11 11:20:20 UTC (rev 896)
@@ -1,24 +0,0 @@
-# CHNOSZ/demo/protein_Cp.R
-# Compare calculated and experimental Cp of proteins
-# 20250510 demo extracted from anintro.Rmd
-
-library(CHNOSZ)
-
-# Experimental values from Privalov and Makhatadze (1990)
-PM90 <- read.csv(system.file("extdata/misc/PM90.csv", package = "CHNOSZ"))
-plength <- protein.length(colnames(PM90)[2:5])
-Cp_expt <- t(t(PM90[, 2:5]) / plength)
-matplot(PM90[, 1], Cp_expt, type = "p", pch = 19,
- xlab = axis.label("T"), ylab = axis.label("Cp0"), ylim = c(110, 280))
-for(i in 1:4) {
- pname <- colnames(Cp_expt)[i]
- aq <- subcrt(pname, "aq", T = seq(0, 150))$out[[1]]
- cr <- subcrt(pname, "cr", T = seq(0, 150))$out[[1]]
- lines(aq$T, aq$Cp / plength[i], col = i)
- lines(cr$T, cr$Cp / plength[i], col = i, lty = 2)
-}
-legend("right", legend = colnames(Cp_expt),
- col = 1:4, pch = 19, lty = 1, bty = "n", cex = 0.9)
-legend("bottomright", legend = c("experimental", "calculated (aq)",
- "calculated (cr)"), lty = c(NA, 1, 2), pch = c(19, NA, NA), bty = "n")
-title("Calculated and experimental Cp of proteins", font.main = 1)
Modified: pkg/CHNOSZ/inst/NEWS.Rd
===================================================================
--- pkg/CHNOSZ/inst/NEWS.Rd 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/inst/NEWS.Rd 2025-05-11 11:20:20 UTC (rev 896)
@@ -15,7 +15,7 @@
\newcommand{\Cp}{\ifelse{latex}{\eqn{C_P}}{\ifelse{html}{\out{<I>C<sub>P</sub></I>}}{Cp}}}
\newcommand{\DG0}{\ifelse{latex}{\eqn{{\Delta}G^{\circ}}}{\ifelse{html}{\out{Δ<I>G</I>°}}{ΔG°}}}
-\section{Changes in CHNOSZ version 2.1.0-67 (2025-05-10)}{
+\section{Changes in CHNOSZ version 2.1.0-68 (2025-05-11)}{
\subsection{OBIGT DEFAULT DATA}{
\itemize{
@@ -117,10 +117,10 @@
\file{CHNOSZ_equilibrium.Rmd} in the [JMDplots
package](https://github.com/jedick/JMDplots).
- \item Add \file{demo/protein_A.R}, \file{demo/protein_Cp.R}, and
- \file{demo/rubisco_Zc.R}, extracted from old version of
+ \item \file{demo/rubisco_Zc.R} extracted from old version of
\file{anintro.Rmd}.
+
}
}
@@ -193,7 +193,7 @@
\itemize{
\item Move \code{read.fasta()}, \code{count.aa()}, and \code{aasum()} to
- canprot package with different names.
+ the canprot package with different names.
\item Remove \code{seq2aa()}.
@@ -200,6 +200,9 @@
\item Remove \code{add.alpha()}. Instead, \code{grDevices::adjustcolor()}
is now used in \code{stack_mosaic()} to add transparency.
+ \item Move \code{protein.equil()} to JMDplots; the associated demo has
+ been moved to the \file{CHNOSZ_equilibrium.Rmd} vignette in that package.
+
}
}
Modified: pkg/CHNOSZ/inst/TODO
===================================================================
--- pkg/CHNOSZ/inst/TODO 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/inst/TODO 2025-05-11 11:20:20 UTC (rev 896)
@@ -6,8 +6,6 @@
- diagram(): produce error or warning for NA values of affinity
-- Move protein.equil() example (previously in equilibrium.Rnw) to JMDplots.
-
- Add elements from Robie and Hemingway (1995).
- Add check to mosaic() that 'stable' values have the right dimensions.
Modified: pkg/CHNOSZ/inst/tinytest/test-protein.info.R
===================================================================
--- pkg/CHNOSZ/inst/tinytest/test-protein.info.R 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/inst/tinytest/test-protein.info.R 2025-05-11 11:20:20 UTC (rev 896)
@@ -10,36 +10,3 @@
info <- "pinfo() returns NA if no organism is given"
expect_equal(pinfo(c("LYSC_CHICK", "MYGPHYCA")), c(6, NA), info = info)
-info <- "protein.equil() reports values consistent with Dick and Shock (2011)"
-protein <- pinfo(c("CSG_METVO", "CSG_METJA"))
-# To reproduce the calculations in the paper, use superseded properties of [Met], [Gly], and [UPBB]
-mod.OBIGT("[Met]", G = -35245, H = -59310, S = 40.38, E_units = "cal")
-mod.OBIGT("[Gly]", G = -6075, H = -5570, S = 17.31, E_units = "cal")
-mod.OBIGT("[UPBB]", G = -21436, H = -45220, S = 1.62, E_units = "cal")
-basis("CHNOS+")
-suppressMessages(swap.basis("O2", "H2"))
-pequil <- protein.equil(protein, loga.protein=-3)
-
-# Before 20230630:
-# With R = 8.314445 (= 1.9872 * 4.184) in util.units(), ionize.aa(pinfo(protein)) gives:
-# 20 18
-#[1,] -56.06509 -55.87025
-# Astar/RT in the paragraph following Eq. 23, p. 6 of DS11
-# (truncated because of rounding)
-# expect_true(any(grepl(c("0\\.435.*1\\.36"), pequil)), info = info)
-
-# After 20230630:
-# With R = 8.314463 in util.units(), ionize.aa(pinfo(protein)) gives:
-# 20 18
-#[1,] -56.06511 -55.87027
-# The differences propagate up, so the test was changed on 20230630:
-expect_true(any(grepl(c("0\\.437.*1\\.36"), pequil)), info = info)
-
-# log10 activities of the proteins in the left-hand column of the same page
-expect_true(any(grepl(c("-3\\.256.*-2\\.834"), pequil)), info = info)
-
-# Reference
-
-# Dick, J. M. and Shock, E. L. (2011) Calculation of the relative chemical stabilities of proteins
-# as a function of temperature and redox chemistry in a hot spring.
-# PLoS ONE 6, e22782. https://doi.org/10.1371/journal.pone.0022782
Modified: pkg/CHNOSZ/man/examples.Rd
===================================================================
--- pkg/CHNOSZ/man/examples.Rd 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/man/examples.Rd 2025-05-11 11:20:20 UTC (rev 896)
@@ -13,13 +13,13 @@
\usage{
examples(save.png = FALSE)
- demos(which = c("references", "protein.equil", "affinity", "NaCl",
+ demos(which = c("references", "affinity", "NaCl",
"density", "ORP", "ionize", "buffer", "protbuff",
"glycinate", "mosaic", "copper", "arsenic", "solubility", "gold",
"contour", "sphalerite", "minsol", "Shh", "saturation",
"adenine", "DEW", "lambda", "potassium", "TCA", "aluminum", "AD",
"comproportionation", "Pourbaix", "E_coli", "yttrium", "rank.affinity",
- "uranyl", "sum_S", "MgATP", "protein_Cp", "protein_A", "rubisco_Zc"),
+ "uranyl", "sum_S", "MgATP", "rubisco_Zc"),
save.png = FALSE)
}
@@ -33,7 +33,6 @@
\describe{
\item{references}{Cross-check the references in refs.csv with the thermodynamic database}
- \item{protein.equil}{Chemical activities of two proteins in metastable equilibrium (Dick and Shock, 2011)}
\item{affinity}{Affinities of metabolic reactions and amino acid synthesis (Amend and Shock, 1998, 2001)}
\item{NaCl}{Equilibrium constant for aqueous NaCl dissociation (Shock et al., 1992)}
\item{density}{Density of \H2O, inverted from IAPWS-95 equations (\code{\link{rho.IAPWS95}})}
@@ -69,8 +68,6 @@
\item{uranyl}{Total (carbonate|sulfate)-pH diagrams for uranyl species (Migdisov et al., 2024)}
\item{sum_S}{Summed molality of S species and solubility contours for iron and gold (Skirrow and Walshe, 2002)}
\item{MgATP}{Speciation of ATP with H+ and Mg+2 (Alberty, 2003)}
- \item{protein_Cp}{Compare calculated and experimental Cp of proteins (Privalov and Makhatadze, 1990)}
- \item{protein_A}{Affinity of ionization of proteins}
\item{rubisco_Zc}{Zc of Rubisco vs optimal growth temperature}
}
@@ -143,8 +140,6 @@
Pourbaix, M. (1974) \emph{Atlas of Electrochemical Equilibria in Aqueous Solutions}, NACE, Houston, TX and CEBELCOR, Brussels. \url{https://www.worldcat.org/oclc/563921897}
-Privalov, P. L. and Makhatadze, G. I. (1990) Heat capacity of proteins. II. Partial molar heat capacity of the unfolded polypeptide chain of proteins: Protein unfolding effects. \emph{J. Mol. Biol.} \bold{213}, 385--391. \doi{10.1016/S0022-2836(05)80198-6}
-
Schulte, M. D. and Shock, E. L. (1995) Thermodynamics of Strecker synthesis in hydrothermal systems. \emph{Orig. Life Evol. Biosph.} \bold{25}, 161--173. \doi{10.1007/BF01581580}
Shock, E. L. and Koretsky, C. M. (1995) Metal-organic complexes in geochemical processes: Estimation of standard partial molal thermodynamic properties of aqueous complexes between metal cations and monovalent organic acid ligands at high pressures and temperatures. \emph{Geochim. Cosmochim. Acta} \bold{59}, 1497--1532. \doi{10.1016/0016-7037(95)00058-8}
Modified: pkg/CHNOSZ/man/ionize.aa.Rd
===================================================================
--- pkg/CHNOSZ/man/ionize.aa.Rd 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/man/ionize.aa.Rd 2025-05-11 11:20:20 UTC (rev 896)
@@ -58,8 +58,50 @@
legend = c("PM90 experiment", "MP90 groups",
"DLH06 groups no ion", "DLH06 groups ionized"))
title("Heat capacity of unfolded LYSC_CHICK")
+
+## Compare calculated and experimental Cp of proteins
+## Using affinity() interface to ionization calculations
+# Experimental values from Privalov and Makhatadze (1990)
+PM90 <- read.csv(system.file("extdata/misc/PM90.csv", package = "CHNOSZ"))
+plength <- protein.length(colnames(PM90)[2:5])
+Cp_expt <- t(t(PM90[, 2:5]) / plength)
+matplot(PM90[, 1], Cp_expt, type = "p", pch = 19,
+ xlab = axis.label("T"), ylab = axis.label("Cp0"), ylim = c(110, 280))
+for(i in 1:4) {
+ pname <- colnames(Cp_expt)[i]
+ aq <- subcrt(pname, "aq", T = seq(0, 150))$out[[1]]
+ cr <- subcrt(pname, "cr", T = seq(0, 150))$out[[1]]
+ lines(aq$T, aq$Cp / plength[i], col = i)
+ lines(cr$T, cr$Cp / plength[i], col = i, lty = 2)
}
+legend("right", legend = colnames(Cp_expt),
+ col = 1:4, pch = 19, lty = 1, bty = "n", cex = 0.9)
+legend("bottomright", legend = c("experimental", "calculated (aq)",
+ "calculated (cr)"), lty = c(NA, 1, 2), pch = c(19, NA, NA), bty = "n")
+title("Calculated and experimental Cp of proteins", font.main = 1)
+## Calculate affinity of ionization of proteins
+ip <- pinfo(c("CYC_BOVIN", "LYSC_CHICK", "MYG_PHYCA", "RNAS1_BOVIN"))
+basis("CHNOS+")
+a_ion <- affinity(pH = c(0, 14), iprotein = ip)
+basis("CHNOS")
+a_nonion <- affinity(iprotein = ip)
+plot(c(0, 14), c(50, 300), xlab = "pH", ylab = quote(italic(A/2.303*RT)), type = "n")
+for(i in 1:4) {
+ A_ion <- as.numeric(a_ion$values[[i]])
+ A_nonion <- as.numeric(a_nonion$values[[i]])
+ lines(a_ion$vals[[1]], A_ion - A_nonion, col=i)
+}
+legend("topright", legend = a_ion$species$name,
+ col = 1:4, lty = 1, bty = "n", cex = 0.9)
+title("Affinity of ionization of proteins", font.main = 1)
+
+# NOTE: The affinity is always positive, showing the strong energetic drive for
+# ionization of proteins in aqueous solution. The degrees of ionization of
+# amino and carboxyl groups increase at low and high pH, respectively, giving
+# rise to the U-shaped lines.
+}
+
\references{
Dick, J. M., LaRowe, D. E. and Helgeson, H. C. (2006) Temperature, pressure, and electrochemical constraints on protein speciation: Group additivity calculation of the standard molal thermodynamic properties of ionized unfolded proteins. \emph{Biogeosciences} \bold{3}, 311--336. \doi{10.5194/bg-3-311-2006}
Modified: pkg/CHNOSZ/man/protein.info.Rd
===================================================================
--- pkg/CHNOSZ/man/protein.info.Rd 2025-05-10 11:02:43 UTC (rev 895)
+++ pkg/CHNOSZ/man/protein.info.Rd 2025-05-11 11:20:20 UTC (rev 896)
@@ -6,7 +6,6 @@
\alias{protein.formula}
\alias{protein.OBIGT}
\alias{protein.basis}
-\alias{protein.equil}
\title{Summaries of thermodynamic properties of proteins}
\description{
@@ -19,7 +18,6 @@
protein.formula(protein, organism = NULL, residue = FALSE)
protein.OBIGT(protein, organism = NULL, state=thermo()$opt$state)
protein.basis(protein, T = 25, normalize = FALSE)
- protein.equil(protein, T=25, loga.protein = 0, digits = 4)
}
\arguments{
@@ -30,8 +28,6 @@
\item{normalize}{logical, return per-residue values (those of the proteins divided by their lengths)?}
\item{state}{character, physical state}
\item{T}{numeric, temperature in \degC}
- \item{loga.protein}{numeric, decimal logarithms of reference activities of proteins}
- \item{digits}{integer, number of significant digits (see \code{\link{signif}})}
}
\details{
@@ -62,12 +58,15 @@
\code{protein.basis} calculates the numbers of the basis species (i.e. opposite of the coefficients in the formation reactions) that can be combined to form the composition of each of the proteins.
The basis species must be present in \code{thermo()$basis}, and if \samp{H+} is among the basis species, the ionization states of the proteins are included.
-The ionization state of the protein is calculated at the pH defined in \code{thermo()$basis} and at the temperature specified by the \code{T} argument.
+The ionization state of the protein is calculated using \code{\link{ionize.aa}} at the pH defined in \code{thermo()$basis} and at the temperature specified by the \code{T} argument.
If \code{normalize} is TRUE, the coefficients on the basis species are divided by the lengths of the proteins.
- \code{protein.equil} produces a series of messages showing step-by-step a calculation of the chemical activities of proteins in metastable equilibrium. For the first protein, it shows the standard Gibbs energies of the reaction to form the nonionized protein from the basis species and of the ionization reaction of the protein (if \samp{H+} is in the basis), then the standard Gibbs energy/RT of the reaction to form the (possibly ionized) protein per residue. The per-residue values of \samp{logQstar} and \samp{Astar/RT} are also shown for the first protein. Equilibrium calculations are then performed, only if more than one protein is specified. This calculation applies the Boltzmann distribution to the calculation of the equilibrium degrees of formation of the residue equivalents of the proteins, then converts them to activities of proteins taking account of \code{loga.protein} and protein length. If the \code{protein} argument is numeric (indicating rownumbers in \code{thermo()$protein}), the values of \samp{Astar/RT} are compared with the output of \code{\link{affinity}}, and those of the equilibrium degrees of formation of the residues and the chemical activities of the proteins with the output of \code{\link{diagram}}. If the values in any of these tests are are not \code{\link{all.equal}} an error is produced indicating a bug.
}
+\seealso{
+ \code{\link{ionize.aa}}
+}
+
\examples{\dontshow{reset()}
# Search by name in thermo()$protein
# These are the same: ip1 == ip2
More information about the CHNOSZ-commits
mailing list