1

我是 R 新手，我试图让曲线适合这个散布数据，给我一个高斯曲线。我真的很感激任何帮助。数据：

library(tidyverse)
MK20 <- tribble(~X.Intensity,    ~Average,
             0.400,  0.0000000,
             0.463,  0.0000000,
             0.536,  0.000000,
             0.621,  0.0000000,
             0.719,  0.0000000,
             0.833,  0.0000000,
             0.965,  0.0000000,
             1.120,  0.0000000,
             1.290,  0.0000000,
             1.500,  0.0000000,
             1.740,  0.0000000,
             2.010,  0.0000000,
             2.330,  0.0000000,
             2.700,  0.0000000,
             3.120,  0.0000000,
             3.620,  0.0000000,
             4.190,  0.0000000,
             4.850,  0.0000000,
             5.610,  0.0000000,
             6.500,  0.0000000,
             7.530,  0.0000000,
             8.720,  0.0000000,
             10.100,  0.0000000,
             11.700,  0.0000000,
             13.500,  0.0000000,
             15.700,  0.0000000,
             18.200,  0.0000000,
             21.000,  0.0000000,
             24.400,  0.0000000,
             28.200,  0.0000000,
             32.700,  0.0000000,
             37.800,  0.0000000,
             43.800,  0.7023333,
             50.700,  3.3700000,
             58.800,  7.3933333,
             68.100, 11.4666667,
             78.800, 14.3666667,
             91.300, 15.4000000,
             106.000, 14.5000000,
             122.000, 12.0000000,
             142.000,  8.6433333,
             164.000,  5.2200000,
             190.000,  2.4500000,
             220.000,  0.7580000,
             255.000,  0.1306667,
             295.000,  0.0000000,
             342.000,  0.0000000,
             396.000,  0.0000000,
             459.000,  0.0000000,
             531.000,  0.0000000,
             615.000,  0.0000000,
             712.000,  0.0000000,
             825.000,  0.0000000,
             955.000,  0.0000000,
             1110.000,  0.0000000,
             1280.000,  0.0000000,
             1480.000,  0.0000000,
             1720.000,  0.0000000,
             1990.000,  0.0000000,
             2300.000,  0.0000000,
             2670.000,  0.0000000,
             3090.000,  0.0000000,
             3580.000,  0.0000000,
             4150.000,  0.0000000,
             4800.000,  0.0000000,
             5560.000,  0.0000000,
             6440.000,  0.0000000,
             7460.000,  0.0000000,
             8630.000,  0.0000000)

我用来绘制的代码是：

plot(log10(MK20$X.Intensity), MK20$Average, col=1, xlim=c(-0.5,4), ylim=c(0,20), xlab="Log(Average diameter)", ylab="Intensity", xaxt='n')

我正在使用 minor.tick.axis 函数在对数 x 轴上添加次要刻度。我想向该数据添加一条高斯曲线（最适合）。我试图type='l'在绘图上添加一个，但曲线并不平滑，我不想要一条必须触及每个数据点但最适合的曲线。

很抱歉，如果解决方案很简单，但我无法弄清楚。

4

2 回答 2

0

在这种情况下，我们不能使用通常的fitdistr方法来拟合正态分布，因为我们没有原始数据。看起来“平均”列是某种类型的密度估计。如果它是 pdf，那么它应该集成到 1 但它没有。

f <- approxfun(x = log10(MK20$X.Intensity), y= MK20$Average)
integrate(f, lower = log10(0.4), upper = log10(8630))

#6.142134 with absolute error < 0.00043

因此，我们可以通过将其缩小约 6.14 来将其转换为 pdf，然后尝试找到与该 pdf 匹配的均值和标准差。

这是一个简单的高斯拟合的第一次尝试。首先，我选择了平均值 2（通过查看密度最大的位置）、k = 6.14（积分值）的比例因子，然后使用 sd 玩，直到有一个合理的拟合。

m=2
s=0.15
k=6.14
x_seq = seq(1,3,length.out = 100)
df <- tibble(x_seq = x_seq, dens = dnorm(x_seq, m, s))


MK20 %>% 
  mutate(log_intensity = log10(X.Intensity)) %>% 
  ggplot(aes(log_intensity, Average/k)) +
  geom_point() +
  geom_line(data = df, aes(x_seq, dens))

在此处输入图像描述

接下来，我使用 optimx 通过最小化拟合和数据之间的平方和来拟合 3 个参数（k = 比例因子，m = 平均值，s = 标准偏差）。

目标函数（拟合和数据之间差异的平方和）

f <- function(x) {
  k = x[1]
  m = x[2]
  s = x[3]
  MK20 %>% 
  mutate(log_intensity = log10(X.Intensity)) %>%
  mutate(fit = dnorm(log_intensity, m, s)) %>% 
  summarise(sum((fit - Average/k)^2)) %>% pull
}

使用 optimx 查找参数（最小平方和） 参数的初始值取自眼睛拟合。

library(optimx)    
optimx(par = c(6.14, 2, 0.15), fn = f )

#k = 6.294696 m = 1.971488 s= 0.1583936

让我们用拟合的参数重新绘制

# points for a gaussian
x_seq = seq(1,3,length.out = 100) 
df <- tibble(x_seq = x_seq, dens = dnorm(x_seq, m, s))


MK20 %>% 
  mutate(log_intensity = log10(X.Intensity)) %>% 
  ggplot(aes(log_intensity, Average/k)) +
  geom_point() +
  geom_line(data = df, aes(x_seq, dens))

在此处输入图像描述

于 2019-12-10T11:05:45.797 回答
0

接触每个点的曲线肯定会最适合您的数据。:)

除此之外，您可以尝试包含平滑曲线，例如

plot(log10(MK20$X.Intensity), MK20$Average, col=1, xlim=c(-0.5,4), ylim=c(0,20), 
     xlab="Log(Average diameter)", ylab="Intensity", xaxt='n', type='n')
lines(lowess(MK20$Average ~ log10(MK20$X.Intensity), f=0.3))

您可以f=在（0 和 1）之间改变参数以更改平滑级别。

这是 f=0.3 时的输出。

在此处输入图像描述

于 2019-12-10T10:22:41.207 回答