【问题标题】:Normalizing histogram bins in gnuplot在gnuplot中标准化直方图箱
【发布时间】:2011-08-12 19:49:57
【问题描述】:

我正在尝试绘制一个直方图,其 bin 由 bin 中的元素数量标准化。

我正在使用以下

binwidth=5
bin(x,width)=width*floor(x/width) + binwidth/2.0
plot 'file' using (bin($2, binwidth)):($4) smooth freq with boxes

得到一个基本的直方图,但我希望每个 bin 的值除以 bin 的大小。我如何在 gnuplot 中解决这个问题,或者在必要时使用外部工具?

【问题讨论】:

    标签: gnuplot histogram normalize bins


    【解决方案1】:

    简单

    plot 'file' using (bin($2, binwidth)):($4/$4) smooth freq with boxes
    

    【讨论】:

    • 您能否在您的答案周围添加一些上下文。它将帮助提问者和读者。
    • 该说什么?这是对shivknight问题的直接回答。通过使用($4/$4)(1) 而不仅仅是($4),可以获得元素的计数而不是它们的总和。
    【解决方案2】:

    计算文件中数据点数量的另一种方法是使用系统命令。如果您要绘制多个文件,并且事先不知道点数,这将非常有用。我用过:

    countpoints(file) = system( sprintf("grep -v '^#' %s| wc -l", file) )
    file1count = countpoints (file1)
    file2count = countpoints (file2)
    file3count = countpoints (file3)
    ...
    

    countpoints 函数避免计算以“#”开头的行。然后,您将使用已经提到的函数来绘制归一化直方图。

    这是一个完整的例子:

    n=100
    xmin=-50.
    xmax=50.
    binwidth=(xmax-xmin)/n
    
    bin(x,width)=width*floor(x/width)+width/2.0
    countpoints(file) = system( sprintf("grep -v '^#' %s| wc -l", file) )
    
    file1count = countpoints (file1)
    file2count = countpoints (file2)
    file3count = countpoints (file3)
    
    plot file1 using (bin(($1),binwidth)):(1.0/(binwidth*file1count)) smooth freq with boxes,\
         file2 using (bin(($1),binwidth)):(1.0/(binwidth*file2count)) smooth freq with boxes,\
         file3 using (bin(($1),binwidth)):(1.0/(binwidth*file3count)) smooth freq with boxes
    ...
    

    【讨论】:

      【解决方案3】:

      在gnuplot 4.6中,可以通过stats命令统计点数,比plot快。其实你不需要s(x)=((sum=sum+1),0)这样的技巧,而是直接在stats 'out.dat' u 1运行后通过变量STATS_records来统计个数。

      【讨论】:

        【解决方案4】:

        在 gnuplot 4.4 中,函数具有不同的属性,因为它们可以执行多个连续的命令,然后返回一个值(参见 gnuplot tricks)这意味着您实际上可以计算点数 n,在gnuplot 文件,而无需事先知道。此代码针对一个文件“out.dat”运行,该文件包含一列:来自正态分布的 n 个样本的列表:

        binwidth = 0.1
        set boxwidth binwidth
        sum = 0
        
        s(x)          = ((sum=sum+1), 0)
        bin(x, width) = width*floor(x/width) + binwidth/2.0
        
        plot "out.dat" u ($1):(s($1))
        plot "out.dat" u (bin($1, binwidth)):(1.0/(binwidth*sum)) smooth freq w boxes
        

        第一个绘图语句读取数据文件并为每个点增加一次 sum,绘制一个零。

        第二个绘图语句实际上是使用 sum 的值对直方图进行归一化。

        【讨论】:

        • 您可以通过让s(x) 的第二个值成为NaN 并将notitle 添加到第一个plot 命令来进一步改进这一点 - 这样,总和将完全不可见在图中,由于 gnuplot 在绘图时忽略了NaN 值 =)
        【解决方案5】:

        我会这样做,使用以下命令从 R 生成 n=500 个随机高斯变量:

        Rscript -e 'cat(rnorm(500), sep="\\n")' > rnd.dat
        

        我使用与您完全相同的想法来定义标准化直方图,其中 y 定义为 1/(binwidth * n),除了我使用 int 而不是 floor 并且我没有在仓值。简而言之,这是对smooth.dem 演示脚本的快速改编,Janert 的教科书中描述了类似的方法,Gnuplot in ActionChapter 13,第 257 页,免费提供)。您可以将我的示例数据文件替换为 random-points,它位于 Gnuplot 附带的 demo 文件夹中。请注意,我们需要将点数指定为 Gnuplot,因为文件中的记录没有计数功能。

        bw1=0.1
        bw2=0.3
        n=500
        bin(x,width)=width*int(x/width)
        set xrange [-3:3]
        set yrange [0:1]
        tstr(n)=sprintf("Binwidth = %1.1f\n", n) 
        set multiplot layout 1,2
        set boxwidth bw1
        plot 'rnd.dat' using (bin($1,bw1)):(1./(bw1*n)) smooth frequency with boxes t tstr(bw1)
        set boxwidth bw2
        plot 'rnd.dat' using (bin($1,bw2)):(1./(bw2*n)) smooth frequency with boxes t tstr(bw2)
        

        这是结果,有两个 bin 宽度

        此外,这确实是一种粗略的直方图方法,在 R 中很容易获得更详细的解决方案。事实上,问题是如何定义一个好的 bin 宽度,这个问题已经在 stats.stackexchange.com 上讨论过:使用 @ 987654325@ 分箱规则实施起来应该不会太难,但您需要计算四分位间距。

        这是 R 将如何处理相同的数据集,使用默认选项(Sturges 规则,因为在这种特殊情况下,这不会产生影响)和上面使用的等间距 bin。

        使用的R代码如下:

        par(mfrow=c(1,2), las=1)
        hist(rnd, main="Sturges", xlab="", ylab="", prob=TRUE)
        hist(rnd, breaks=seq(-3.5,3.5,by=.1), main="Binwidth = 0.1", 
             xlab="", ylab="", prob=TRUE)
        

        您甚至可以通过检查调用 hist() 时返回的值来了解 R 是如何工作的:

        > str(hist(rnd, plot=FALSE))
        List of 7
         $ breaks     : num [1:14] -3.5 -3 -2.5 -2 -1.5 -1 -0.5 0 0.5 1 ...
         $ counts     : int [1:13] 1 1 12 20 49 79 108 87 71 43 ...
         $ intensities: num [1:13] 0.004 0.004 0.048 0.08 0.196 0.316 0.432 0.348 0.284 0.172 ...
         $ density    : num [1:13] 0.004 0.004 0.048 0.08 0.196 0.316 0.432 0.348 0.284 0.172 ...
         $ mids       : num [1:13] -3.25 -2.75 -2.25 -1.75 -1.25 -0.75 -0.25 0.25 0.75 1.25 ...
         $ xname      : chr "rnd"
         $ equidist   : logi TRUE
         - attr(*, "class")= chr "histogram"
        

        这就是说,如果您愿意,您可以使用 R 结果通过 Gnuplot 处理您的数据(尽管我建议直接使用 R :-)。

        【讨论】:

          猜你喜欢
          • 2017-08-28
          • 2020-10-26
          • 1970-01-01
          • 1970-01-01
          • 2022-01-20
          • 2016-02-25
          • 2016-01-13
          • 2020-12-04
          • 1970-01-01
          相关资源
          最近更新 更多