前言
之前写过一篇绘制博客园积分与排名趋势图的文章——《查看博客园积分与排名趋势图的工具 》，使用那篇文章介绍的工具，可以通过趋势图直观的看出排名前进的走势。但是如果想看看自己积分达到多少才能进入前多少名次，就无能为力了。如果我们能够根据历史数据，拟合出一条预测曲线，然后根据这条曲线就可以预测多少积分进多少排名啦！想想就很激动呐~
积分-排名曲线在开始拟合数据之前，我们先看下现在的趋势图：
它的横轴是时间轴，两条纵轴分别是积分与排名。虽然积分总是增长、排名大体是下降趋势，但是这个趋势我试了几种方法都不太好拟合，感觉没什么规律，特别受发表文章时机的影响：如果每天发一篇，这个曲线肯定增长的快；如果一年也不发，增长的肯定慢。我们的预测其实和时间没有什么关系，主要自变量是积分，因变量是排名。那么是否可以做一张图，横轴是积分，纵轴是排名呢？于是就有了下面的改进：
这下清爽多了，可以看到，虽然发表文章时积分会快速增长从而在横轴留下一些空白，但是总的来看积分与排名相对趋势是维持在一条曲线上的。想要绘制这样一条曲线，gnuplot 脚本改动并不大：
代码语言：javascript代码运行次数：0运行复制 1 #! /usr/bin/gnuplot 2 set terminal png size 1080,720 #建立空白图片 3 set title usr.": score (".y1range.") rank (".y2range.")" #注明曲线图标题 4 set output "fit.png" #设置文件名 5 set key left 6 set grid 7 8 set xlabel "score" 9 set ylabel "rank"10 11 #plot "score.txt" using 1:2 with lp pt 13 title "score" axis x1y1, "score.txt" using 1:3 with lp pt 3 title "rank" axis x1y212 plot "score.txt" using 2:3 with lp pt 13 title "score-rank"13 14 quit #退出软件就是将 line 11 换成 line 12，再去掉一些冗余代码就可以了。
数据拟合有了历史数据和正确的映射关系，就可以进行数据拟合了。数据拟合最重要的是找到拟合函数，第一眼看到上面那条曲线我想到的就是二次函数，可以用抛物线的一段来进行拟合。此外随着积分的无限增大，排名最终会无限贴进于 1，这个特性让我想到了很多无限贴进某个值的函数，例如倒数函数无限贴进于 0，对数函数的反函数无限贴近于 0，还有反正切函数也是如此。最近刚刚高考结束，相信大家对高中数学已经忘了个七七八八，还好我这里准备了几张图：
从左到右依次是倒数、对数、反正切(atan)，对于后两者，还需要加工一下才能满足我们的要求，具体罗列如下：
代码语言：javascript代码运行次数：0运行复制f1(x)=a*x^2+b*x+cf2(x)=f/x+gf3(x)=j*atan(x)+kf4(x)=m*log(x)+n
分别就是它们四个啦，为了之后区分参数方便，使用了不同的参数名称。有的人可能会问了，这个对数是无限趋近于无穷大的，怎么能和上面的曲线拟合在一起呢？反正切也存在类似的问题，先别着急，让我们看下拟合结果就知道了。
嘿，这样居然也行，四条曲线竟然都能拟合上。gnuplot 脚本改动并不大：
代码语言：javascript代码运行次数：0运行复制 1 #! /usr/bin/gnuplot 2 set terminal png size 1080,720 #建立空白图片 3 set title usr.": score (".y1range.") rank (".y2range.")" #注明曲线图标题 4 set output "fit.png" #设置文件名 5 set key left 6 set grid 7 8 set xlabel "score" 9 set ylabel "rank"10 11 y1(x)=a*x**2+b*x+c12 fit y1(x) "score.txt" using 2:3 via a,b,c13 14 y2(x)=f/x+g15 fit y2(x) "score.txt" using 2:3 via f,g16 17 y3(x)=j*atan(x)+k18 fit y3(x) "score.txt" using 2:3 via j,k19 20 y4(x)=m*log(x)+n21 fit y4(x) "score.txt" using 2:3 via m,n22 23 plot "score.txt" using 2:3 with lp pt 13 title "score-rank", \24 y1(x) with l lw 2 lt 2 title "f(x)=ax^2+bx+c", \25 y2(x) with l lw 3 lt 3 title "f(x)=a/x+b", \26 y3(x) with l lw 1 lt 4 title "f(x)=a*atan(x)+b", \27 y4(x) with l lw 2 lt 5 title "f(x)=a*log(x)+b"28 29 quit #退出软件line 11-21 定义了各个拟合函数，line 24-27 增加了拟合曲线的绘制。如果能将拟合后的函数参数标识出来，就更好了，其实也不难，因为 a/b/c/f/g/j/k/m/n 这些参数在 gnuplot 脚本中就可以直接访问，只需要在图例显示处增加一些代码就可以了：
代码语言：javascript代码运行次数：0运行复制plot "score.txt" using 2:3 with lp pt 13 title "score-rank", \ y1(x) with l lw 6 lt 2 title sprintf("f1(x)=%.8fx^2%+fx%+.0f",a,b,c), \ y2(x) with l lw 5 lt 3 title sprintf("f2(x)=%.2f/x%+.0f",f,g), \ y3(x) with l lw 2 lt 4 title sprintf("f3(x)=%.2fatan(x)%+.0f",j,k), \ y4(x) with l lw 2 lt 5 title sprintf("f4(x)=%.2flog(x)%+.0f",m,n)最终效果如下：
看到对数函数的参数，可以解答之前的疑问了，m=-39186.57 (数据预测
有了三个拟合函数，就可以对数据进行预测了，一开始雄心勃勃，打算预测一下自己 40 W 分时的排名 (有点不自量力哈)，预测值通过 label 形式输出在图形上，就像这样：
结果相去甚远，首先恭喜二次函数 (f1) 得到了 2800 W 名的好成绩~ 话说我刚加入博客园也才 14 W 名，简直离谱；倒数函数得到了 4500 名的结果，比现在 (5900) 也好不了多少，这 38 W 分看来白涨了； 对数函数最牛，直接干成负数了，祝贺我进入了另一个维度~ 看来还是我好高骛远了，手里拿着 2 W 分的数据，却想预测 40 W 分的排名，难为 gnuplot 了。调整了一下目标，先算个 4 W 分的吧：
现在稍微靠谱一些了，总的来说，二次函数偏差太大，倒数函数偏大，对数函数偏小。其实想想也能明白，二次函数抛物线嘛，迟早要转回来的，偏大也可以理解。最后上 gnuplot 脚本：
代码语言：javascript代码运行次数：0运行复制 1 #! /usr/bin/gnuplot 2 set terminal png size 1080,720 #建立空白图片 3 set title usr.": score (".y1range.") rank (".y2range.")" #注明曲线图标题 4 set output "fit.png" #设置文件名 5 set key left 6 set grid 7 8 set xlabel "score" 9 set ylabel "rank"10 11 y1(x)=a*x**2+b*x+c12 fit y1(x) "score.txt" using 2:3 via a,b,c13 14 y2(x)=f/x+g15 fit y2(x) "score.txt" using 2:3 via f,g16 17 y3(x)=m*log(x)+n18 fit y3(x) "score.txt" using 2:3 via m,n19 20 xval=4000021 y1val=a*xval**2+b*xval+c22 y2val=f/xval+g23 y3val=m*log(xval)+n24 25 set label 1 sprintf("f1(%.0f)=%.0f",xval,y1val) at graph 0.6,0.7 left26 set label 2 sprintf("f2(%.0f)=%.0f",xval,y2val) at graph 0.6,0.65 left27 set label 3 sprintf("f3(%.0f)=%.0f",xval,y3val) at graph 0.6,0.6 left28 29 plot "score.txt" using 2:3 with lp pt 13 title "score-rank", \30 y1(x) with l lw 4 lt 2 title sprintf("f1(x)=%.8fx^2%+fx%+.0f",a,b,c), \31 y2(x) with l lw 3 lt 3 title sprintf("f2(x)=%.2f/x%+.0f",f,g), \32 y3(x) with l lw 2 lt rgb "red" title sprintf("f3(x)=%.2flog(x)%+.0f",m,n)33 34 35 quit #退出软件主要增加的内容就是中间几行啦，line 20-23 通过现有参数计算 x 为 40000 时的各函数 y 值，line 25-27 将计算好的预测值显示在 label 中。其实函数已经定义好了，如果能直接通过 f1 (40000) / f2 (40000) / f3 (40000) 得到结果就更好了，但是没有在 gnuplot 手册中找到这种语法，不得己自己再写一遍，有懂行的同学不吝赐教哈。最后因为我们的预测值都是整数，所以打印出来的数据也没有保留小数位，通过 sprintf 自动四舍五入了。
绘制预测曲线上面的代码可以预测某个点的数据，但是还是有点呆板，需要手动指定预测值，如果将预测值设置为当前分数的两倍，就能自动预测啦。将得到的预测值写入一个数据文件，随着时间积累，形成一条预测曲线绘制出来，再和实际数据做对比，预测效果岂不一目了然？
输出预测值将 gnuplot 脚本中计算得到的预测值写入一个文件，这个事情看起来简单做起来难，难就难在我找了半天，没有找到可以从脚本直接输出信息到 console 或重定向到文件的方法。echo 这种命令在 gnuplot 脚本中是不存在的，于是这里绕了一个大圈——在脚本执行完成后，通过分拆 fit.log 中的拟合日志提取函数的各个参数 (a/b/c/f/g/m/n)，再构建函数计算预测值，最后写入数据文件——哪位高手如果知道如何在 gnuplot 脚本中直接输出信息的话，不吝赐教哈，就可以把这个大弯路省掉了。
现在来看这个蹩脚的弯路，在开始参数提取前，先熟悉一下 gnuplot 的拟合日志：
代码语言：javascript代码运行次数：0运行复制*******************************************************************************Mon Jul 5 14:22:31 2021FIT: data read from "score.txt" using 2:3 format = x:z #datapoints = 365 residuals are weighted equally (unit weight)function used for fitting: y1(x)y1(x)=a*x**2+b*x+cfitted parameters initialized with current variable valuesiter chisq delta/lim lambda a b c 0 1.2882895936e+19 0.00e+00 1.08e+08 1.000000e+00 1.000000e+00 1.000000e+00 11 6.0109804999e+08 -9.02e-01 1.08e-03 1.991252e-04 -8.363164e+00 1.465348e+05After 11 iterations the fit converged.final sum of squares of residuals : 6.01098e+08rel. change during last iteration : -9.01929e-06degrees of freedom (FIT_NDF) : 362rms of residuals (FIT_STDFIT) = sqrt(WSSR/ndf) : 1288.6variance of residuals (reduced chisquare) = WSSR/ndf : 1.66049e+06Final set of parameters Asymptotic Standard Error======================= ==========================a = 0.000199125 +/- 3.57e-06 (1.793%)b = -8.36316 +/- 0.08565 (1.024%)c = 146535 +/- 456.8 (0.3117%)correlation matrix of the fit parameters: a b c a 1.000 b -0.986 1.000 c 0.921 -0.969 1.000 *******************************************************************************Mon Jul 5 14:22:31 2021FIT: data read from "score.txt" using 2:3 format = x:z #datapoints = 365 residuals are weighted equally (unit weight)function used for fitting: y2(x)y2(x)=f/x+gfitted parameters initialized with current variable valuesiter chisq delta/lim lambda f g 0 2.5365572521e+12 0.00e+00 7.07e-01 1.000000e+00 1.000000e+00 7 2.5241195880e+09 -4.95e-08 7.07e-08 3.478261e+08 4.432964e+04After 7 iterations the fit converged.final sum of squares of residuals : 2.52412e+09rel. change during last iteration : -4.94761e-13degrees of freedom (FIT_NDF) : 363rms of residuals (FIT_STDFIT) = sqrt(WSSR/ndf) : 2636.95variance of residuals (reduced chisquare) = WSSR/ndf : 6.9535e+06Final set of parameters Asymptotic Standard Error======================= ==========================f = 3.47826e+08 +/- 2.776e+06 (0.7982%)g = 44329.6 +/- 327.3 (0.7383%)correlation matrix of the fit parameters: f g f 1.000 g -0.907 1.000 *******************************************************************************Mon Jul 5 14:22:31 2021FIT: data read from "score.txt" using 2:3 format = x:z #datapoints = 365 residuals are weighted equally (unit weight)function used for fitting: y3(x)y3(x)=m*log(x)+nfitted parameters initialized with current variable valuesiter chisq delta/lim lambda m n 0 2.5360128666e+12 0.00e+00 6.58e+00 1.000000e+00 1.000000e+00 5 2.4184679129e+08 -5.46e-07 6.58e-05 -3.918657e+04 4.438066e+05After 5 iterations the fit converged.final sum of squares of residuals : 2.41847e+08rel. change during last iteration : -5.46492e-12degrees of freedom (FIT_NDF) : 363rms of residuals (FIT_STDFIT) = sqrt(WSSR/ndf) : 816.238variance of residuals (reduced chisquare) = WSSR/ndf : 666245Final set of parameters Asymptotic Standard Error======================= ==========================m = -39186.6 +/- 95.82 (0.2445%)n = 443807 +/- 886.9 (0.1998%)correlation matrix of the fit parameters: m n m 1.000 n -0.999 1.000
它将拟合的迭代、最终的参数、误差率清楚的罗列了出来，这里我们需要提取的只是各个参数最终的值。好在这些行都比较有特点，基本遵循 "参数名 = 数据 ***" 的格式，于是可以用 grep 先过滤一把：
代码语言：javascript代码运行次数：0运行复制sed -n '/[abcfgmn] *=.*/p' fit.log
经过这样处理后，就只剩下这样的内容了：
代码语言：javascript代码运行次数：0运行复制a = 0.000199125 +/- 3.57e-06 (1.793%)b = -8.36316 +/- 0.08565 (1.024%)c = 146535 +/- 456.8 (0.3117%)f = 3.47826e+08 +/- 2.776e+06 (0.7982%)g = 44329.6 +/- 327.3 (0.7383%)m = -39186.6 +/- 95.82 (0.2445%)n = 443807 +/- 886.9 (0.1998%)
再使用 awk 提取我们关心的前三列：
代码语言：javascript代码运行次数：0运行复制sed -n '/[abcfgmn] *=.*/p' fit.log | awk '{print $1,$2,$3+0}'注意第三列使用 "$3+0" 的 trick 来保证提取的是浮点数据：
代码语言：javascript代码运行次数：0运行复制a = 0.000199125b = -8.36316c = 146535f = 347826000g = 44329.6m = -39186.6n = 443807
最后将它们传递给 eval 实现“求值”，即当前 shell 中自动获得了 a-n 7 个变量和它们的初值：
代码语言：javascript代码运行次数：0运行复制eval $(sed -n '/[abcfgmn] *=.*/p' fit.log | awk '{print $1,$2,$3+0}' | sed 's/ //g')注意 shell 中的赋值语句是不能有空格的，所以需要使用 sed 过滤一下空格。有了函数的各个参数，就可以利用函数计算预测值了，这里我使用了 awk，因为它对浮点数运算支持的比较好：
代码语言：javascript代码运行次数：0运行复制awk -v a=$a -v b=$b -v c=$c -v f=$f -v g=$g -v m=$m -v n=$n -v xval=$xval 'BEGIN { print "y1="int(a*xval*xval+b*xval+c+0.5); print "y2="int(f/xval+g+0.5); print "y3="int(m*log(xval)+n+0.5) }'这里首先利用 awk 的 -v 选项将 shell 脚本中的变量传递到 awk 中，然后在 awk 中根据三个函数分别计算了三个预测值。xval 是事先定义好的，一般就是当前 x 值的 2 倍。上面的脚本输出如下：
y1=130609y2=53025y3=28561
这里的四舍五入使用了 +0.5 的笨办法，最终结果和 gnuplot 计算的完全一致。有了这个就可以继续利用 eval 搞事情了，下面上完整的提取、计算、写入脚本：
代码语言：javascript代码运行次数：0运行复制 1 #! /bin/sh 2 usr=$(cat user.txt) 3 firstline=$(cat ./score.txt | head -1) 4 y1min=$(echo $firstline | awk '{ print $2 }') 5 y2max=$(echo $firstline | awk '{ print $3 }') 6 lastline=$(cat ./score.txt | tail -1) 7 y1max=$(echo $lastline | awk '{ print $2 }') 8 y2min=$(echo $lastline | awk '{ print $3 }') 9 echo "y1 range [$y1min,$y1max], y2 range [$y2max,$y2min]"10 gnuplot -e "usr='$usr'" -e "y1range='$(($y1max-$y1min))'" -e "y2range='$(($y2max-$y2min))'" ./draw.plt11 12 # clear log to get generated values13 rm fit.log 2>/dev/null14 gnuplot -e "usr='$usr'" -e "y1max='$y1max'" -e "y2min='$y2min'" ./fit.plt15 16 if [ $# -gt 0 ]; then 17 # don't know how to print out parameter (a,b,c,f,g,m,n) from gnuplot script, so here extract them from fit.log (that's why we clear fit.log before calling fit.plt)18 eval $(sed -n '/[abcfgmn] *=.*/p' fit.log | awk '{print $1,$2,$3+0}' | sed 's/ //g')19 # after that line, a/b/c/f/g/m/n takes effect, now calculate predicating values20 21 # predicate x*2, round result to integer22 xval=$(($y1max*2))23 eval $(awk -v a=$a -v b=$b -v c=$c -v f=$f -v g=$g -v m=$m -v n=$n -v xval=$xval 'BEGIN { print "y1="int(a*xval*xval+b*xval+c+0.5); print "y2="int(f/xval+g+0.5); print "y3="int(m*log(xval)+n+0.5) }')24 25 # dump results to data files26 echo "$xval $y1" >> predicate_binomial.data27 echo "$xval $y2" >> predicate_reciprocal.data28 echo "$xval $y3" >> predicate_logarithm.data29 fi30 31 # for centos32 type eog > /dev/null 2>&133 if [ $? -eq 0 ]; then 34 eog draw.png &35 eog fit.png &36 exit 037 fi38 39 # for mac40 type open > /dev/null 2>&141 if [ $? -eq 0 ]; then 42 open draw.png &43 open fit.png &44 exit 045 fi46 47 # for windows msys248 type mspaint > /dev/null 2>&149 if [ $? -eq 0 ]; then 50 mspaint draw.png &51 mspaint fit.png &52 exit 053 fi54 55 exit 1重点就是 line 16-29 啦，这段代码只有在脚本参数个数大于 0 时才生效，也就是说你平时做一个 ./plot.sh 这种动作查看趋势图时，是不会修改数据的，只有当定时任务每日执行 score.sh 才会调用带参数的 plot.sh 去更新预测值。关于 score.sh 的内容，可以参数我之前写的那篇文章。
预测值经过计算并提取到 shell 脚本后，分别存储在了三个 data 文件中，文件名说明了他们使用的拟合函数。这里也可以将多个拟合函数的预测值放在一个文件，毕竟他们的 x 轴数据都是一样的嘛，没有这样做的原因主要是考虑到后期可能加入新的拟合函数来进行预测，独立存储的话互不干扰，加入删除都比较方便，利于扩展。
历史补算以现在 20000 分得到的预测值是在 40000 分左右，想要验证预测值准不准，还需要我的博客涨到 40000 分才能拉在一起做对比，以我现在这种速度，不知道要等到猴年马月，然而我是这种没有耐心的人吗？确实是！于是我打算将 5000-20000 这段历史数据的预测值也给它补全了，这样生成的预测值就涵盖了 10000-40000 的范围，正好能和我 10000-20000 的区间形成交集，从而对比着看。
说起来复杂，做起来简单，我只需将历史数据一分为二，一部分是 4000-5000 的极小规模的历史数据；一部分是 5000-20000 的补算历史数据。使用 plot.sh 作用于第一部分数据，生成预测值，然后从第二部分数据头部取出一条记录添加到第一部分数据末尾，再调用 plot.sh 生成一条预测数据……周而复始，直到第二部分数据消耗完毕。下面就是用于补算的脚本：
代码语言：javascript代码运行次数：0运行复制1 #! /bin/sh 2 # data should be parted into score.txt (with some basic historical data) 3 # and more.txt (with more recent data), 4 # and then we will repair predicating data by moving more.txt data 5 # into score.txt line by line and calling plot.sh without call png starter... 6 7 if [ ! -f "more.txt" ]; then 8 echo "you should split score.txt to score.txt & more.txt first before run this scripts..." 9 exit 110 fi11 12 while read line 13 do14 echo "repair line $line"15 echo "$line" >> score.txt16 ./plot.sh "update_predicating_data"17 done < more.txt18 19 echo "repair done, now submit predicat_*.data !"20 rm more.txt
当然了，在运行这个脚本之前，你需要将 score.txt 一切为二，前一部分还叫 score.txt，第二部分需要命名为 more.txt。当脚本运行完毕后，more.txt 中的数据自动合入 score.txt，同时产生三个包含历史预测值的 predicate_xxx.data 文件。而且为了防止在补算历史过程中自动打开大量的趋势图文件，需要在 plot.sh 中做如下修改：
代码语言：javascript代码运行次数：0运行复制 1 if [ $# -gt 0 ]; then 2 # don't know how to print out parameter (a,b,c,f,g,m,n) from gnuplot script, so here extract them from fit.log (that's why we clear fit.log before calling fit.plt) 3 eval $(sed -n '/[abcfgmn] *=.*/p' fit.log | awk '{print $1,$2,$3+0}' | sed 's/ //g') 4 # after that line, a/b/c/f/g/m/n takes effect, now calculate predicating values 5 6 # predicate x*2, round result to integer 7 xval=$(($y1max*2)) 8 eval $(awk -v a=$a -v b=$b -v c=$c -v f=$f -v g=$g -v m=$m -v n=$n -v xval=$xval 'BEGIN { print "y1="int(a*xval*xval+b*xval+c+0.5); print "y2="int(f/xval+g+0.5); print "y3="int(m*log(xval)+n+0.5) }') 9 10 # dump results to data files11 echo "$xval $y1" >> predicate_binomial.data12 echo "$xval $y2" >> predicate_reciprocal.data13 echo "$xval $y3" >> predicate_logarithm.data14 else15 # for centos16 type eog > /dev/null 2>&117 if [ $? -eq 0 ]; then 18 eog draw.png &19 eog fit.png &20 exit 021 fi22 23 # for mac24 type open > /dev/null 2>&125 if [ $? -eq 0 ]; then 26 open draw.png &27 open fit.png &28 exit 029 fi30 31 # for windows msys232 type mspaint > /dev/null 2>&133 if [ $? -eq 0 ]; then 34 mspaint draw.png &35 mspaint fit.png &36 exit 037 fi38 fi39 40 exit 1line 14 添加了一行 else，表示只有在脚本参数为 0 时才启动图片自动打开功能。
绘制预测线前面铺垫了这么多，终于可以把预测值绘制出来一睹芳容了：
先撇开预测曲线的风骚走位，重点关注一下 10000-20000 这个区间，可以看到点划线的真实数据和三条预测曲线相差还是蛮大的，只有对数函数在一开始非常贴近真实值，后来也出现偏离。不过总的来说对数最准、倒数整体偏高、二次函数就是来搞笑的——上下横跳没准头。有的人可能奇怪了，这个预测值为什么和拟合曲线差距这么大？原因是预测曲线的每个点的参数都不一样，由之前小一半的历史数据拟合计算得到的，所以不能完美重合拟合函数，可以将预测曲线理解成是一堆拟合函数的末位点集合形成的轨迹 (稍费脑，理解不了就不用理解了)。最后再来欣赏一下全图，无意间拉大横轴范围后，二次函数的抛物线本性果然暴露了，对数和倒数表现倒还可以。过足了预测瘾后，还是让我们限制一下横轴范围，不然真实数据实在是看不清了：
同时调整的还有图例显示方式和预测曲线的颜色，后者可以让你一眼看出拟合函数与预测曲线关系。下面是最终的 gnuplot 脚本：
代码语言：javascript代码运行次数：0运行复制 1 #! /usr/bin/gnuplot 2 set terminal png size 1080,720 #建立空白图片 3 set title usr.": score (".y1max.") rank (".y2min.")" #注明曲线图标题 4 set output "fit.png" #设置文件名 5 set key left reverse Left spacing 1.2 6 set grid 7 8 set xlabel "score" 9 set ylabel "rank"10 # to prevent predicating value pollute our x-axis11 set xrange [y1min-100:y1max+100]12 # no fit log in console (fit.log only)13 set fit quiet14 15 y1(x)=a*x**2+b*x+c16 fit y1(x) "score.txt" using 2:3 via a,b,c17 18 y2(x)=f/x+g19 fit y2(x) "score.txt" using 2:3 via f,g20 21 y3(x)=m*log(x)+n22 fit y3(x) "score.txt" using 2:3 via m,n23 24 xval=y1max*225 y1val=a*xval**2+b*xval+c26 y2val=f/xval+g27 y3val=m*log(xval)+n28 29 set label 1 sprintf("f1(%.0f)=%.0f",xval,y1val) at graph 0.6,0.7 left30 set label 2 sprintf("f2(%.0f)=%.0f",xval,y2val) at graph 0.6,0.65 left31 set label 3 sprintf("f3(%.0f)=%.0f",xval,y3val) at graph 0.6,0.6 left32 33 plot "score.txt" using 2:3 with lp pt 13 title "score-rank", \34 y1(x) with l lw 4 lt 2 title sprintf("f1(x)=%.8fx^2%+fx%+.0f",a,b,c), \35 y2(x) with l lw 3 lt 3 title sprintf("f2(x)=%.2f/x%+.0f",f,g), \36 y3(x) with l lw 2 lt rgb "red" title sprintf("f3(x)=%.2flog(x)%+.0f",m,n), \37 "predicate_binomial.data" using 1:2 with lp pt 12 lt 2 title "f1-pred", \38 "predicate_reciprocal.data" using 1:2 with lp pt 11 lt 3 title "f2-pred", \39 "predicate_logarithm.data" using 1:2 with lp pt 10 lt rgb "red" title "f3-pred"40 41 quit #退出软件主要改动点在 line 37-39，使用之前补算的预测值绘图。
结论最后比较靠谱的预测值就是取对数函数与倒数函数预测值靠近对数函数 1/3 的位置，即 f(x)=(f3(x)-f2(x))/3+f2(x)=2/3*f2(x)+1/3*f3(x)。
这个就是我胡诌的啦~ 真正的结论就是自己的数值分析没学好，需要回炉重造，连一个好的拟合函数也找不到，此刻终于在毕业 N 多年后认识到高等数据和数值分析的重要性…
后记文中代码可以从这里下载：
https://github.com/goodpaperman/cnblogs
这个代码库可以直接 fork 哟~ 在写文章过程中发现了一个很有意思的数学教学网站，可以在 web 端直接绘制你输入的表达试，推荐一波：https://www.shuxuele.com/index.html
参考1. gnuplot图例legend设置
2. awk将字符串转为数字的方法
3. 在命令行中使用gnuplot快速查看数据
4. Gnuplot重定向fit输出
5. gnuplot常用技巧
6. 在gnuplot中，绘制一些分段函数
7. gnuplot使用手册
8. shell脚本，awk实现跳过文件里面的空行。
9. AWK 打印匹配内容之后的指定行
10. UWP 获取博客园积分，并以图表形式呈现变化趋势
