Bio.SeqUtils.IsoelectricPoint module
Calculate isoelectric points of polypeptides using methods of Bjellqvist.
pK values and the methods are taken from:
* Bjellqvist, B.,Hughes, G.J., Pasquali, Ch., Paquet, N., Ravier, F.,
Sanchez, J.-Ch., Frutiger, S. & Hochstrasser, D.F.
The focusing positions of polypeptides in immobilized pH gradients can be
predicted from their amino acid sequences. Electrophoresis 1993, 14,
1023-1031.
* Bjellqvist, B., Basse, B., Olsen, E. and Celis, J.E.
Reference points for comparisons of two-dimensional maps of proteins from
different human cell types defined in a pH scale where isoelectric points
correlate with polypeptide compositions. Electrophoresis 1994, 15, 529-539.
I designed the algorithm according to a note by David L. Tabb, available at: http://fields.scripps.edu/DTASelect/20010710-pI-Algorithm.pdf
- class Bio.SeqUtils.IsoelectricPoint.IsoelectricPoint(protein_sequence, aa_content=None)
Bases:
object
A class for calculating the IEP or charge at given pH of a protein.
- Parameters:
- :protein_sequence: A ``Bio.Seq`` or string object containing a protein
sequence.
- :aa_content: A dictionary with amino acid letters as keys and its
occurrences as integers, e.g.
{"A": 3, "C": 0, ...}
. Default:None
. IfNone
, the dic will be calculated from the given sequence.
Examples
The methods of this class can either be accessed from the class itself or from a
ProtParam.ProteinAnalysis
object (with partially different names):>>> from Bio.SeqUtils.IsoelectricPoint import IsoelectricPoint as IP >>> protein = IP("INGAR") >>> print(f"IEP of peptide {protein.sequence} is {protein.pi():.2f}") IEP of peptide INGAR is 9.75 >>> print(f"Its charge at pH 7 is {protein.charge_at_pH(7.0):.2f}") Its charge at pH 7 is 0.76
>>> from Bio.SeqUtils.ProtParam import ProteinAnalysis as PA >>> protein = PA("PETER") >>> print(f"IEP of {protein.sequence}: {protein.isoelectric_point():.2f}") IEP of PETER: 4.53 >>> print(f"Charge at pH 4.53: {protein.charge_at_pH(4.53):.2f}") Charge at pH 4.53: 0.00
Methods
:charge_at_pH(pH): Calculates the charge of the protein for a given pH
:pi(): Calculates the isoelectric point
- __init__(protein_sequence, aa_content=None)
Initialize the class.
- charge_at_pH(pH)
Calculate the charge of a protein at given pH.
- pi(pH=7.775, min_=4.05, max_=12)
Calculate and return the isoelectric point as float.
This is a recursive function that uses bisection method. Wiki on bisection: https://en.wikipedia.org/wiki/Bisection_method
- Arguments:
pH: the pH at which the current charge of the protein is computed. This pH lies at the centre of the interval (mean of min_ and max_).
min_: the minimum of the interval. Initial value defaults to 4.05, which is below the theoretical minimum, when the protein is composed exclusively of aspartate.
max_: the maximum of the the interval. Initial value defaults to 12, which is above the theoretical maximum, when the protein is composed exclusively of arginine.