Histograms and multimodal distribution detection, fixes #429

samples/BenchmarkDotNet.Samples/Intro/IntroMultimodal.cs

@@ -0,0 +1,24 @@
+using System;
+using System.Threading;
+using BenchmarkDotNet.Attributes;
+using BenchmarkDotNet.Attributes.Columns;
+using BenchmarkDotNet.Attributes.Jobs;
+using BenchmarkDotNet.Engines;
+
+namespace BenchmarkDotNet.Samples.Intro
+{
+    [MValueColumn]
+    [SimpleJob(RunStrategy.Throughput, 1, 0, -1, 1, "MainJob")]
+    public class IntroMultimodal
+    {
+        private readonly Random rnd = new Random(42);
+
+        private void Multimodal(int n) 
+            => Thread.Sleep((rnd.Next(n) + 1) * 100);
+
+        [Benchmark] public void Unimodal() => Multimodal(1);
+        [Benchmark] public void Bimodal() => Multimodal(2);
+        [Benchmark] public void Trimodal() => Multimodal(3);
+        [Benchmark] public void Quadrimodal() => Multimodal(4);
+    }
+}

src/BenchmarkDotNet.Core/Analysers/MultimodalDistributionAnalyzer.cs

@@ -0,0 +1,37 @@
+using System.Collections.Generic;
+using System.Linq;
+using BenchmarkDotNet.Engines;
+using BenchmarkDotNet.Extensions;
+using BenchmarkDotNet.Jobs;
+using BenchmarkDotNet.Mathematics;
+using BenchmarkDotNet.Reports;
+using JetBrains.Annotations;
+
+namespace BenchmarkDotNet.Analysers
+{
+    public class MultimodalDistributionAnalyzer : AnalyserBase
+    {
+        public override string Id => "MultimodalDistribution";
+        public static readonly IAnalyser Default = new MultimodalDistributionAnalyzer();
+
+        private MultimodalDistributionAnalyzer() { }
+
+        [NotNull]
+        private Conclusion Create([NotNull] string kind, double mValue, [CanBeNull] BenchmarkReport report)
+            => CreateWarning($"It seems that the distribution {kind} (mValue = {mValue})", report);
+
+        public override IEnumerable<Conclusion> AnalyseReport(BenchmarkReport report, Summary summary)
+        {
+            var statistics = report.ResultStatistics;
+            if (statistics == null || statistics.N < 15)
+                yield break;
+            double mValue = MathHelper.CalculateMValue(statistics);
+            if (mValue > 4.2)
+                yield return Create("is multimodal", mValue, report);
+            else if (mValue > 3.2)
+                yield return Create("is bimodal", mValue, report);
+            else if (mValue > 2.8)
+                yield return Create("can have several modes", mValue, report);
+        }
+    }
+}

src/BenchmarkDotNet.Core/Attributes/Columns/IterationsColumnAttribute.cs

@@ -0,0 +1,11 @@
+using BenchmarkDotNet.Columns;
+using JetBrains.Annotations;
+
+namespace BenchmarkDotNet.Attributes.Columns
+{
+    [PublicAPI]
+    public class IterationsColumnAttribute : ColumnConfigBaseAttribute
+    {
+        public IterationsColumnAttribute() : base(StatisticColumn.Iterations) { }
+    }
+}

src/BenchmarkDotNet.Core/Attributes/Columns/MValueColumnAttribute.cs

@@ -0,0 +1,17 @@
+using BenchmarkDotNet.Columns;
+using JetBrains.Annotations;
+
+namespace BenchmarkDotNet.Attributes.Columns
+{
+    /// <summary>
+    /// Prints mvalue.
+    /// See http://www.brendangregg.com/FrequencyTrails/modes.html
+    /// </summary>
+    [PublicAPI]
+    public class MValueColumnAttribute: ColumnConfigBaseAttribute
+    {
+        public MValueColumnAttribute() : base(StatisticColumn.MValue)
+        {
+        }
+    }
+}

src/BenchmarkDotNet.Core/Columns/StatisticColumn.cs

@@ -17,7 +17,7 @@ public enum Priority
             Percentiles,
             Additional
         }
-
+        
         public static readonly IColumn Mean = new StatisticColumn("Mean", "Arithmetic mean of all measurements",
             s => s.Mean, Priority.Main);
 
@@ -53,6 +53,15 @@ public enum Priority
         public static readonly IColumn Kurtosis = new StatisticColumn("Kurtosis", "Measure of the tailedness ( fourth standardized moment)",
             s => s.Kurtosis, Priority.Additional, UnitType.Dimensionless);
 
+        /// <summary>
+        /// See http://www.brendangregg.com/FrequencyTrails/modes.html
+        /// </summary>
+        public static readonly IColumn MValue = new StatisticColumn("MValue", "Modal value, see http://www.brendangregg.com/FrequencyTrails/modes.html",
+            MathHelper.CalculateMValue, Priority.Additional, UnitType.Dimensionless);
+
+        public static readonly IColumn Iterations = new StatisticColumn("Iterations", "Number of target iterations",
+            s => s.N, Priority.Additional, UnitType.Dimensionless);
+
         public static readonly IColumn P0 = CreatePercentileColumn(0, s => s.Percentiles.P0);
         public static readonly IColumn P25 = CreatePercentileColumn(25, s => s.Percentiles.P25);
         public static readonly IColumn P50 = CreatePercentileColumn(50, s => s.Percentiles.P50);

src/BenchmarkDotNet.Core/Configs/DefaultConfig.cs

@@ -45,6 +45,7 @@ public IEnumerable<IAnalyser> GetAnalysers()
             yield return OutliersAnalyser.Default;
             yield return MinIterationTimeAnalyser.Default;
             yield return IterationSetupCleanupAnalyser.Default;
+            yield return MultimodalDistributionAnalyzer.Default;
         }
 
         public IEnumerable<IValidator> GetValidators()

src/BenchmarkDotNet.Core/Exporters/PlainExporter.cs

@@ -27,7 +27,7 @@ public override void ExportToLog(Summary summary, ILogger logger)
                 {
                     logger.WriteLine();
                     logger.WriteLineHeader($"* Statistics for {mode}");
-                    logger.WriteLineStatistic(runs.Where(it => it.IterationMode == mode).GetStatistics().ToTimeStr());
+                    logger.WriteLineStatistic(runs.Where(it => it.IterationMode == mode).GetStatistics().ToTimeStr(calcHistogram: true));
                 }
             }
         }

src/BenchmarkDotNet.Core/Extensions/StatisticsExtensions.cs

@@ -1,6 +1,7 @@
 using System.Text;
 using BenchmarkDotNet.Horology;
 using BenchmarkDotNet.Mathematics;
+using BenchmarkDotNet.Mathematics.Histograms;
 
 namespace BenchmarkDotNet.Extensions
 {
@@ -24,7 +25,7 @@ public static string ToStr(this Statistics s)
             return builder.ToString();
         }
 
-        public static string ToTimeStr(this Statistics s, TimeUnit unit = null)
+        public static string ToTimeStr(this Statistics s, TimeUnit unit = null, bool calcHistogram = false)
         {
             if (s == null)
                 return NullSummaryMessage;
@@ -34,12 +35,20 @@ public static string ToTimeStr(this Statistics s, TimeUnit unit = null)
             string errorPercent = (s.StandardError / s.Mean * 100).ToStr("0.00");
             var ci = s.ConfidenceInterval;
             string ciMarginPercent = (ci.Margin / s.Mean * 100).ToStr("0.00");
+            double mValue = MathHelper.CalculateMValue(s);
             builder.AppendLine($"Mean = {s.Mean.ToTimeStr(unit)}, StdErr = {s.StandardError.ToTimeStr(unit)} ({errorPercent}%); N = {s.N}, StdDev = {s.StandardDeviation.ToTimeStr(unit)}");
             builder.AppendLine($"Min = {s.Min.ToTimeStr(unit)}, Q1 = {s.Q1.ToTimeStr(unit)}, Median = {s.Median.ToTimeStr(unit)}, Q3 = {s.Q3.ToTimeStr(unit)}, Max = {s.Max.ToTimeStr(unit)}");
             builder.AppendLine($"IQR = {s.InterquartileRange.ToTimeStr(unit)}, LowerFence = {s.LowerFence.ToTimeStr(unit)}, UpperFence = {s.UpperFence.ToTimeStr(unit)}");
             builder.AppendLine($"ConfidenceInterval = {s.ConfidenceInterval.ToTimeStr(unit)}, Margin = {ci.Margin.ToTimeStr(unit)} ({ciMarginPercent}% of Mean)");
-            builder.AppendLine($"Skewness = {s.Skewness.ToStr()}, Kurtosis = {s.Kurtosis.ToStr()}");
-            return builder.ToString();
+            builder.AppendLine($"Skewness = {s.Skewness.ToStr()}, Kurtosis = {s.Kurtosis.ToStr()}, MValue = {mValue.ToStr()}");
+            if (calcHistogram)
+            {
+                var histogram = HistogramBuilder.Adaptive.Build(s);
+                builder.AppendLine("-------------------- Histogram --------------------");
+                builder.AppendLine(histogram.ToTimeStr());
+                builder.AppendLine("---------------------------------------------------");
+            }
+            return builder.ToString().Trim();
         }
     }
 }

src/BenchmarkDotNet.Core/Mathematics/Histograms/AdaptiveHistogramBuilder.cs

@@ -0,0 +1,162 @@
+using System;
+using System.Collections.Generic;
+using System.Linq;
+using BenchmarkDotNet.Extensions;
+using JetBrains.Annotations;
+
+namespace BenchmarkDotNet.Mathematics.Histograms
+{
+    public class AdaptiveHistogramBuilder : IHistogramBuilder
+    {
+        [PublicAPI, Pure]
+        public Histogram Build(Statistics s, BinSizeRule? rule = null)
+        {
+            double binSize = s.GetOptimalBinSize(rule);
+            if (Math.Abs(binSize) < 1e-9)
+                binSize = 1;
+            return BuildWithFixedBinSize(s.GetValues(), binSize);
+        }
+
+        // TODO: Optimize
+        [PublicAPI, Pure]
+        public Histogram BuildWithFixedBinSize(IEnumerable<double> values, double binSize)
+        {
+            const double eps = 1e-9;
+            const double margin = 0.1;
+            const double adaptiveFactor = 0.02;
+
+            if (binSize < eps)
+                throw new ArgumentException($"binSize ({binSize.ToStr()}) should be a positive number", nameof(binSize));
+
+            var list = values.ToList();
+            if (list.IsEmpty())
+                throw new ArgumentException("Values should be non-empty", nameof(values));
+
+            list.Sort();
+
+            var points = new List<double> { NiceFloor(list.Min() - binSize / 2), NiceCeiling(list.Max() + binSize / 2) };
+            int processedPointCount = 0;
+            while (true)
+            {
+                int pointIndex = -1;
+                for (int i = processedPointCount; i < points.Count - 1; i++)
+                {
+                    double adaptiveBinSize = (points[i] + points[i + 1]) / 2.0 * adaptiveFactor;
+                    double maxSize = Math.Max(binSize * (1.0 + 2 * margin), adaptiveBinSize);
+                    if (points[i + 1] - points[i] > maxSize)
+                    {
+                        pointIndex = i;
+                        break;
+                    }
+                }
+
+                if (pointIndex == -1)
+                    break;
+
+                double lower = points[pointIndex];
+                double upper = points[pointIndex + 1];
+
+                int bestIndex1 = -1;
+                int bestIndex2 = -1;
+                int bestCount = -1;
+                double bestDist = double.MaxValue;
+
+                bool Inside(double x) => x > lower - eps && x < upper - eps;
+
+                for (int i = 0; i < list.Count; i++)
+                    if (Inside(list[i]))
+                    {
+                        int j = i;
+                        while (j < list.Count && Inside(list[j]) && list[j] - list[i] < binSize)
+                            j++;
+                        int count = j - i;
+                        var dist = list[j - 1] - list[i];
+                        if (count > bestCount || count == bestCount && dist < bestDist)
+                        {
+                            bestCount = count;
+                            bestIndex1 = i;
+                            bestIndex2 = j - 1;
+                            bestDist = dist;
+                        }
+                    }
+
+                if (bestIndex1 != -1)
+                {
+                    double center = (list[bestIndex1] + list[bestIndex2]) / 2.0;
+                    double adaptiveBinSize = Math.Max(binSize, center * adaptiveFactor);
+                    double left = center - adaptiveBinSize / 2;
+                    double right = Math.Min(center + adaptiveBinSize / 2, upper);
+
+                    if (left > lower + binSize * margin)
+                        points.Insert(pointIndex + 1, NiceFloor(left));
+                    else if (right < upper - binSize * margin)
+                    {
+                        points.Insert(pointIndex + 1, NiceFloor(right));
+                        processedPointCount++;
+                    }
+                    else
+                        processedPointCount++;
+                }
+                else
+                {
+                    points.Insert(pointIndex + 1, NiceFloor(lower + binSize));
+                    processedPointCount++;
+                }
+            }
+
+            var bins = new List<HistogramBin>(points.Count - 1);
+            int counter = 0;
+            for (int i = 0; i < points.Count - 1; i++)
+            {
+                var bin = new List<double>();
+                double lower = points[i];
+                double upper = points[i + 1];
+
+                while (counter < list.Count && (list[counter] < upper || i == points.Count - 1))
+                    bin.Add(list[counter++]);
+
+                bins.Add(new HistogramBin(lower, upper, bin.ToArray()));
+            }
+
+            // Trim
+            while (bins.Any() && bins.First().IsEmpty)
+                bins.RemoveAt(0);
+            while (bins.Any() && bins.Last().IsEmpty)
+                bins.RemoveAt(bins.Count - 1);
+
+            // Join small bins to neighbors
+            counter = 0;
+            double lastValue = 0;
+            while (counter < bins.Count)
+            {
+                if (bins[counter].HasAny)
+                    lastValue = Math.Max(lastValue, bins[counter].Values.Last());
+                double adaptiveThreshold = Math.Max(binSize / 2, lastValue * adaptiveFactor);
+                if (bins[counter].Gap < adaptiveThreshold)
+                {
+                    double leftGap = counter > 0 ? bins[counter - 1].Gap : double.MaxValue;
+                    double rightGap = counter < bins.Count - 1 ? bins[counter + 1].Gap : double.MaxValue;
+                    if (leftGap < rightGap && counter > 0)
+                    {
+                        bins[counter - 1] = HistogramBin.Union(bins[counter - 1], bins[counter]);
+                        bins.RemoveAt(counter);
+                    }
+                    else if (counter < bins.Count - 1)
+                    {
+                        bins[counter] = HistogramBin.Union(bins[counter], bins[counter + 1]);
+                        bins.RemoveAt(counter + 1);
+                    }
+                    else
+                        counter++;
+                }
+                else
+                    counter++;
+            }
+
+            return new Histogram(binSize, bins.ToArray());
+        }
+
+        private static double NiceFloor(double value) => Math.Floor(value * 1000) / 1000;
+        private static double NiceCeiling(double value) => Math.Ceiling(value * 1000) / 1000;
+    }
+}

src/BenchmarkDotNet.Core/Mathematics/Histograms/BinSizeRule.cs

@@ -0,0 +1,17 @@
+namespace BenchmarkDotNet.Mathematics.Histograms
+{
+    public enum BinSizeRule
+    {
+        FreedmanDiaconis,
+
+        Scott,
+
+        Scott2,
+
+        SquareRoot,
+
+        Sturges,
+
+        Rice,
+    }
+}

src/BenchmarkDotNet.Core/Mathematics/Histograms/Histogram.cs

@@ -0,0 +1,29 @@
+using System.Linq;
+using JetBrains.Annotations;
+
+namespace BenchmarkDotNet.Mathematics.Histograms
+{
+    [PublicAPI]
+    public class Histogram
+    {
+        [PublicAPI]
+        public double BinSize { get; }
+
+        [PublicAPI, NotNull]
+        public HistogramBin[] Bins { get; }
+
+        internal Histogram(double binSize, [NotNull] HistogramBin[] bins)
+        {
+            BinSize = binSize;
+            Bins = bins;
+        }
+
+        // For unit tests
+        [Pure, NotNull]
+        internal static Histogram BuildManual(double binSize, [NotNull] double[][] bins)
+        {
+            var histogramBins = bins.Select(bin => new HistogramBin(bin.Any() ? bin.Min() : 0, bin.Any() ? bin.Max() : 0, bin)).ToArray();
+            return new Histogram(binSize, histogramBins);
+        }
+    }
+}