R统计学(07): 常见数学函数

作者: R语言和Python学堂 | 来源:发表于2019-11-04 22:13 被阅读0次

在统计计算中，你只要学习三类数学规则：分别与幂、指数和对数有关。

对于形如 $y=x^b$ 的表达式，其中x是自变量，指数b为常数，称为幂函数(power function)。
对于形如 $y=e^x$ 的表达式，其中e是自然常数(=2.718282)，指数x为自变量，称为指数函数(exponential function)。
$y=log(x)$ 是指数函数 $y=e^x$ 的反函数(inverse function)，称为对数函数(logarithmic function)。注意：这里的log()是以自然常数e为底的，也即log(x)等同于ln(x)。在R中，log()就是以e为底的对数函数，R中不存在ln()函数。看个例子：

> e <- exp(1)      ### e为自然常数
> e
[1] 2.718282

> log(e)   ### log()函数默认以e为底
[1] 1

> ln(e)    ### R中不存在ln()函数
Error in ln(e) : could not find function "ln"

> log(100, base=10)   ### 通过设置base参数，可以改变对数底的大小，这里为10
[1] 2

更多有关对数的计算可参考这篇教程：R语言初级教程(04): 算术运算

记住一些关于极限的数学结论也是很有用的。我们想知道当x趋于无穷大时y会怎样，以及当x趋于0时y将怎样。下面是一些最重要的结论:

任何数的0次幂都是1： $x^0=1$
1的任何次幂都是1： $1^x=1$
无穷大加1还是无穷大： $\infty + 1 = \infty$
无穷大的倒数为0： $\frac{1}{\infty} = 0$
任何大于1的数的无穷次幂都为无穷大： $1.2^{\infty}=\infty$
若 $0 \leq x < 1$ ，x的无穷次幂为0： $0.999^\infty = 0$
负幂是倒数： $x^{-b} = \frac{1}{x^b}$
自然对数的底为e，值为2.71828，因此： $e^{\infty}=\infty$
分数次幂是开根： $x^{1/3} = \sqrt[3]{x}$
最后一个，但是最有用的一个： $e^{-\infty}=\frac{1}{e^\infty}=\frac{1}{\infty}=0$

下面介绍各种函数：

1. 指数函数(exponential function)和对数函数(Logarithmic function)

指数函数一般表示为： $y=ae^{bx}$ 。

指数函数的反函数为对数函数，对数函数一般表示为： $y=aln(bx)$ ，这里对数的底为自然常数e(=2.71828)。

指数函数和对数函数都是平滑函数，要在R中绘制平滑函数，需要生成一系列在min(x)和max(x)之间的等差数列（元素个数大于100或更多）：

x <- seq(0, 10, 0.01)

在R中，指数函数为exp()，自然对数函数为log()。令a=b=1，指数函数和对数函数的曲线为：

par(mfrow = c(1, 2))
y <- exp(x)
plot(y ~ x, type="l", main="指数函数")
y <- log(x)
plot(y ~ x, type="l", main="对数函数")

2. 三角函数(Trigonometric function)

这里的x(以弧度为单位)在0到2π范围内变化，它的余弦函数(cosine, 底边/斜边)、正弦函数(sine, 垂直边/斜边)和正切函数(tangent, 垂直边/底边)如下图。回想一下，整个圆是2π弧度，所以1弧度 = 360/2π = 57.29578°。

par(mfrow=c(2, 2))
x <- seq(0, 2*pi, 2*pi/10000)
y1 <- cos(x)
y2 <- sin(x)
y3 <- tan(x)
plot(y1~x, type="l", main="余弦")
plot(y2~x, type="l", main="正弦")
plot(y3~x, type="l", ylim=c(-10,10), main="正切")
range(y3)
[1] -1.591549e+03  1.633124e+16

x的正切值是不连续的。在x向右趋于π/2和3π/2时，它的值趋于正无穷大；在x向左趋于π/2和3π/2时，它的值趋于负无穷大。因此，限制y轴上绘制的数值范围(这里从–10到+ 10)可以更好地描绘tan函数的形状。请注意，在ylim定义的范围内，R在x = π/2和x = 3π/2用直线连接正无穷大和负无穷大这两点。

3. 幂率(Power laws)

幂率是一个非常重要的双参数数学函数族，其形式为：

$y=ax^b$

根据幂值b的不同，这种关系可以有五种形式。在b = 0的平凡情况下，函数是y = a(一条水平直线)。下面是四个更有趣的图形：

par(mfrow=c(2, 2))
x <- seq(0, 2, 0.01)
y <- x^0.5    ## a = 1; b = 0.5
plot(x, y, type="l", main="0<b<1")
y <- x    ## a = 1; b = 1
plot(x, y, type="l", main="b=1")
y <- x^2    ## a = 1; b = 2
plot(x, y, type="l", main="b>1")
y <- 1/x    ## a = 1; b = -1
plot(x, y, type="l", main="b<0")

这些函数在许多学科中都很有用。由于该函数经双对数变换(log–log transformation)后变为线性的了，因此参数a和b很容易从数据中估计出来：

$log(y) = log(ax^b) = log(a) + b log(x)$

因此，在对数-对数轴上，截距为log(a)，斜率为b。

泰勒幂定律(Taylor’s power law)是生态昆虫学的一个重要经验关系。它描述了样本的平均值和方差之间的关系。在基本统计模型中，方差被假设为常数(即方差不取决于平均值)。然而，在野外数据中，泰勒发现种群数目(Y)的方差随着平均值µ以幂率的形式增加（即 $var(Y)=a\mu^b$ ），在双对数轴上，大多数系统的数据都落在通过原点斜率为1的直线(泊松分布数据显示的模式，其方差等于平均值，详见R统计学(05): 泊松分布)上方和通过原点斜率为2的直线下方。泰勒幂定律指出，对于特定的系统：

log(variance)是log(mean)的线性函数
log(variance)对log(mean)的回归斜率大于1且小于2，斜率变化范围很小
双对数回归的参数值a和b是系统的基本特征

4. 多项式函数(Polynomial functions)

多项式函数是自变量x出现几次的函数，每次上升到不同的幂。它们有助于描述带有驼峰(hump)、拐点(inflection)或局部极大值(local maximum)的曲线，如：

绘制上图的代码为：

x <- seq(0, 10, 0.1)
y1 <- 2+5*x-0.2*x^2
y2 <- 2+5*x-0.4*x^2
y3 <- 2+4*x-0.6*x^2+0.04*x^3
y4 <- 2+4*x+2*x^2-0.6*x^3+0.04*x^4
par(mfrow=c(2, 2))
plot(x, y1, type="l", ylab="y", main="decelerating")
plot(x, y2, type="l", ylab="y", main="humped")
plot(x, y3, type="l", ylab="y", main="inflection")
plot(x, y4, type="l", ylab="y", main="local maximum")

左上方图显示了减速(斜率减小)函数，其为二次多项式： $y_1=2+5x-0.2x^2$

增大 $x^2$ 项的负系数会产生一个带有驼峰的曲线，如右上方图所示： $y_2=2+5x-0.4x^2$

三次多项式可以显示拐点，如左下方图所示： $y_3=2+4x-0.6x^2+0.04x^3$

最后，含有4次项的多项式能够产生具有局部极大值的曲线，如右下方图所示： $y_4=2+4x+2x^2-0.6x^3+0.04x^4$

多项式的倒数是一类重要的函数，适用于建立广义线性模型(generalized linear models)：
$y= \frac{1}{a + bx + cx^2 + dx^3 + . . . + zx^n}$

根据多项式的阶数(最高次幂)和参数的正负，可以生成各种形状的函数：

x <- seq(0, 10, 0.1)
par(mfrow=c(2, 2))
y1 <- x/(2+5*x)
y2 <- 1/(x-2+4/x)
y3 <- 1/(x^2-2+4/x)
plot(x, y1, type="l", ylab="y", main="Michaelis-Menten")
plot(x, y2, type="l", ylab="y", main="shallow hump")
plot(x, y3, type="l", ylab="y", main="steep hump")

有两种方法可以参数化米氏-米恩方程(Michaelis–Menten equation)： $y=\frac{ax}{1+bx}$ 和 $y=\frac{x}{c+dx}$

在第一种情况下，y的渐近值是a/b，在第二种情况下是1/d。

5. 伽马函数(Gamma function)

伽马函数是阶乘函数 $t!$ 到正实数的扩展：
$\Gamma(t)=\int_0^\infty {x^{t-1}e^{-x}} \,{\rm d}x$

看起来像这样：

绘制上图的代码为：

par(mfrow=c(1, 1))
t <- seq(0.2, 4, 0.01)
plot(t, gamma(t), type="l")
abline(h=1, lty=2)

注意： $\Gamma(1)=\Gamma(2)=1$ ；对于整数t有 $\Gamma(t+1)=t!$

6. 渐近函数(Asymptotic functions)

最常见的渐近函数为：

$y=\frac{ax}{1+bx}$

这个函数在几乎每门学科都有不同的名称。例如，在生物化学中，它被称为米氏-米恩方程(Michaelis–Menten equation)，表示反应速率是酶浓度的函数；在生态学中，它被称为霍林圆盘方程(Holling’s disc equation)，表示捕食者的摄食率是猎物密度的函数。该曲线通过原点并以递减的幅度上升到一个渐近值(x的增加不会导致y进一步增加)。

另一个常见函数是渐近指数：

$y = a(1 − e^{-bx})$

这也是一个双参数模型，在很多情况下，这两个函数对数据的描述所起的效果差不多。

让我们看看这两个渐近函数的极限行为。首先看渐进指数函数，当 $x=0$ 时，有 $y = a(1 − e^{-b×0}) = a(1 − e^0) = 0$ ，因此此函数经过原点。当 $x = \infty$ 时，有 $y = a(1 − e^{-b×\infty}) = a(1 − e^{-\infty}) = a(1 − 0) = a$ ，这证明这个函数是渐进函数，渐进值为a。

接着来看米氏-米恩方程的极限行为。当 $x=0$ 时， $y=\frac{a × 0}{1 + b × 0}=0$ ，此函数也经过原点。当 $x = \infty$ 时，分子分母同时除以x有 $y=\frac{ax}{1+bx}=\frac{a}{\frac{1}{x}+b}=\frac{a}{\frac{1}{\infty}+b}=\frac{a}{b}$ ，因此米氏-米恩方程也是渐进函数，渐进值为a/b。

7. 渐进函数的参数估计

对于渐进指数函数，无法线性化，因此我们必须求助于非线性最小二乘(non-linear least squares, nls)来为其估计参数，后续我们再介绍非线性最小二乘的使用。而米氏-米恩函数的优点之一是易于线性化，使用倒数变换有

$\frac{1}{y}=\frac{1+bx}{ax}=\frac{1}{ax}+\frac{bx}{ax}=\frac{1}{ax}+\frac{b}{a}$

如果令 $Y=1/y$ ， $X=1/x$ ， $A=1/a$ ，以及 $C=b/a$ ，上式子变为：

$Y=AX+C$

这是线性的：C是截距，A是斜率。所以为了从数据中估计a和b的值，我们将x和y都转换成倒数，绘制1/y对1/x的曲线图，进行线性回归，然后反向转换有：

$a=\frac{1}{A}$

$b=aC$

假设我们知道曲线经过两个点(0.2, 44.44)和(0.6, 70.59)。那么我们怎么算出参数a和b的值呢？首先，我们计算四个倒数，线性化函数的斜率A是1/y的变化除以1/x的变化，即：

> (1/44.44 - 1/70.59)/(1/0.2 - 1/0.6)
[1] 0.002500781

因此 $a = 1/A = 1/0.0025 = 400$ 。现在我们重新整理方程，使用其中一个点(比如x = 0.2, y = 44.44)来得到b的值：

$b=\frac{1}{x}(\frac{ax}{y}-1)=\frac{1}{0.2}(\frac{400 × 0.2}{44.44}-1)=4$

8. S型函数(Sigmoid function)

最简单的S型函数是两参数逻辑(two-parameter logistic)函数，它的取值范围 $0<y<1$ 。它的形式为：

$y=\frac{e^{a+bx}}{1+e^{a+bx}}$

这函数对拟合广义线性模型(generalized linear models)至关重要（后续介绍）。

令 $a=1, b=0.5$ ，函数变为： $y=\frac{e^{1+0.5x}}{1+e^{1+0.5x}}$ ，其曲线为：

x <- seq(-10, 10, 0.01)
y <- exp(1+0.5*x)/(1+exp(1+0.5*x))
plot(x, y, type="l", main="two-parameter logistic")

三参数逻辑(three-parameter logistic)函数允许y在任意尺度上变化：

$y=\frac{a}{1+be^{-cx}}$

其截距为 $a/(1+b)$ ，渐进值为 $a$ 。

令 $a=100, b=90, c=1$ ，函数变为： $y=\frac{100}{1+90e^{-x}}$ ，其曲线为：

x <- seq(0, 10, 0.1)
y <- 100/(1+90*exp(-x))
plot(x, y, type="l", main="three-parameter logistic")

四参数逻辑( four-parameter logistic)函数在x轴的左端和右端都有渐近线：

$y=a+\frac{b-a}{1+e^{c(d-x)}}$

其左渐近值为 $a$ ，右渐近值为 $b$ 。

令 $a=20, b=120, c=0.8, d=3$ ，函数变为： $y=20+\frac{100}{1+e^{0.8×(3-x)}}$ ，其曲线为：

x <- seq(-10, 10, 0.1)
y <- 20+100/(1+exp(0.8*(3-x)))
plot(x, y, ylim=c(0,140), type="l", main="four-parameter logistic")

一种不对称S形曲线——贡佩斯生长模型(Gompertz growth model)常用于人口统计和人寿保险工作，它的形式为：

$y = ae^{be^{cx}}$

函数的形状取决于参数 $b$ 和 $c$ 的符号。对于负的S形曲线， $b$ 为负， $c$ 为正；对于正的S形曲线，参数 $b$ 和 $c$ 都为负。看个例子：

## 负的S形曲线, b=-1, c=0.02
par(mfrow=c(1, 2))
x <- -200:100
y <- 100*exp(-exp(0.02*x))
plot(x, y, type="l", main="negative Gompertz")

## 正的S形曲线, b=-5, c=-0.08
x <- 0:100
y <- 50*exp(-5*exp(-0.08*x))
plot(x, y, type="l", main="positive Gompertz")

9. 双指数模型(Biexponential model)

这是一个有用的四参数非线性函数，它是 $x$ 的两个指数函数的和：

$y = ae^{bx} + ce^{dx}$

各种形状取决于参数 $b$ 、 $c$ 和 $d$ (假设 $a$ 为正)的符号。

左上方图显示 $c$ 为正， $b$ 和 $d$ 为负(这是两条指数衰减曲线的总和，因此快速分解材料先消失，然后才是缓慢分解材料，产生两个不同的相位)；
右上方图显示 $c$ 和 $d$ 为正， $b$ 为负(这产生了不对称的U形曲线)；
左下方图显示 $c$ 为负， $b$ 为正， $d$ 为正(有时会产生一条有驼峰的曲线)；
右下方图显示 $b$ 和 $c$ 为正， $d$ 为负。
当 $b$ 、 $c$ 和 $d$ 都为负(未画图)时，该函数被称为一阶隔室模型(first-order compartment model)，用于描述药物在体内的过程，其动力学受三个生理过程影响:代谢、吸收和排泄。

上图的代码为：

x <- seq(0, 10, 0.1)
par(mfrow=c(2, 2))
#1
a <- 10
b <- -0.8
c <- 10
d <- -0.05
y <- a*exp(b*x)+c*exp(d*x)
plot(x, y, main="+ - + -", type="l")
#2
a <- 10
b <- -0.8
c <- 10
d <- 0.05
y <- a*exp(b*x)+c*exp(d*x)
plot(x, y, main="+ - + +", type="l")
#3
a <- 200
b <- 0.2
c <- -1
d <- 0.7
y <- a*exp(b*x)+c*exp(d*x)
plot(x, y, main="+ + - +", type="l")
#4
a <- 200
b <- 0.05
c <- 300
d <- -0.5
y <- a*exp(b*x)+c*exp(d*x)
plot(x, y, main="+ + + -", type="l")

10. 自变量和因变量的转换

我们已经看到了可以使用变换来线性化自变量和因变量之间的关系：

对于指数关系，log(y) ~ x
对于幂函数， log(y) ~ log(x)
对于对数关系， exp(y) ~ x
对于渐近关系， 1/y ~ 1/x

其他转换对于方差稳定很有用：

$\sqrt y$ 稳定计数数据的方差
$arcsin(y)$ 稳定百分比数据的方差

常见数学函数的介绍就到此结束，希望对大家的学习有所帮助，也希望大家多多支持本公众号。

感谢您的阅读！想了解更多有关技巧，请关注我的微信公众号“R语言和Python学堂”，我将定期更新相关文章。

R统计学(07): 常见数学函数

1. 指数函数(exponential function)和对数函数(Logarithmic function)

2. 三角函数(Trigonometric function)

3. 幂率(Power laws)

4. 多项式函数(Polynomial functions)

5. 伽马函数(Gamma function)

6. 渐近函数(Asymptotic functions)

7. 渐进函数的参数估计

8. S型函数(Sigmoid function)

9. 双指数模型(Biexponential model)

10. 自变量和因变量的转换

相关文章

网友评论

延伸阅读

深度阅读

栏目导航

热点阅读

R统计学

R统计学

R与统计分布