def calc_real_underlying_subst_rate(obs_aa_ident_full_protein, rand_ident_region1, rand_ident_region2, fraction_region1_residues=0.3):
    """Calculation of the real, underlying substitution rate in a protein with two regions of differing random identity.

    Estimation of the underlying substitution rate is the first step in calculating a normalisation factor,
    to remove the region1/region2 conservation attributable to the restricted subset of amino acids.

    Please see the associated article manuscript (in preparation) for a full explanation of the method.

    Example calculation:
    region1 = transmembrane(TM)
    region2 = non-transmembrane(nonTM)
    i = obs_aa_ident_full_protein = 0.60  #(i.e. 60% amino acid identity, 40% substitutions)
    fraction_TM_residues = 0.1
    # the observed aa subst rate is 1 - obs_aa_ident_full_protein
    b = 1 - i = 1 - 0.60 = 0.4
    # proportion of seq that is Membranous
    m = fraction_TM_residues = 0.1
    # proportion of seq that is Soluble
    s = 1 - m = 1 - 0.1 = 0.90      # i.e. 90% of the residues are nonTM
    # random identity of TM region
    t = rand_ident_region1 = rand_TM = 0.12
    # random identity of nonTM region
    n = rand_ident_region2 = rand_nonTM = 0.06
    # real underlying aa subst rate for full protein
    # solved from b = mx - mtx + sx - snx
    x = b / ((m*-t) + m - n*s + s) = 0.428

    Parameters
    ----------
    obs_aa_ident_full_protein : float
        Observed amino acid identity of the full pairwise alignment
        "observed" in that it always includes some position that have had a substitution to the same AA (e.g. Leu to Leu), which are unobserved.
        Typically excluding gaps, which are heavily biased by splice variants, etc, and are excluded from AAIMON calculations.
        E.g. 0.6, for 60% amino acid identity.
    rand_TM : float
        Random identity of the TM region.
        e.g. 0.124 for 12.4%
    rand_nonTM : float
        Random identity of the nonTM region.
        e.g. 0.059 for 5.9%
    fraction_TM_residues : float
        Fraction of TM residues in protein sequence, e.g. 0.1 (10% TM, 90% nonTM)
        Used to estimate the real underlying AA substitution rate from the observed AA subst. rate.

    Returns
    -------
    x : float
        Amino acid propensity normalisation factor.
        Real, underlying substitution rate for the full protein.

    Usage
    -----
    # calculate the TM/nonTM conservation ratio for a particular query and homologue
    TM_over_nonTM_cons_ratio = percentage_identity_of_TMD / percentage_identity_of_nonTMD

    # calculate the amino acid propensity normalisation factor.
    obs_aa_ident_full_protein, rand_TM, rand_nonTM, fraction_TM_residues = 0.60, 0.136, 0.055, 0.10
    real_underlying_subst_rate = calc_real_underlying_subst_rate(obs_aa_ident_full_protein, rand_TM, rand_nonTM, fraction_TM_residues)

    aa_prop_norm_factor = calc_aa_prop_norm_factor(real_underlying_subst_rate, rand_TM, rand_nonTM)

    # to normalise, divide the AAIMON by the aa_prop_norm_factor
    TM_over_nonTM_cons_ratio_normalised = TM_over_nonTM_cons_ratio / aa_prop_norm_factor
    """
    # the oBserved aa subst rate is 1 - obs_aa_ident_full_protein
    b = 1 - obs_aa_ident_full_protein
    # proportion of seq that is Membranous
    m = fraction_region1_residues
    # proportion of seq that is Soluble
    s = 1 - m
    # random identity of Tm region
    t = rand_ident_region1
    # random identity of NonTM region
    n = rand_ident_region2
    # real underlying aa subst rate for full protein
    # solved from b = mx - mtx + sx - snx
    x = b / ((m * -t) + m - n * s + s)

    return x