布点算法的原理和实现
在數(shù)據(jù)可視化的過程中，繪制網(wǎng)絡拓撲圖是很重要的，它能清晰呈現(xiàn)一個復雜網(wǎng)絡的結構，節(jié)點的重要性和關系。比如下面幾張圖：
下面這張圖是我的軟件繪制的：
這些都有一個共同的問題，就是如何讓圖繪制的更加美觀？ 復雜的圖結構，已經(jīng)沒法人工布局了。所以計算機自動布局，就成了一個有趣而重要的話題。我們將其稱為布點算法。
基本原理
一個好的布點算法應該能盡量滿足以下四個特點：
- 對稱性：繪制網(wǎng)絡對應將具有相同結構的子圖圍繞繪圖中心平衡布局
- 連接角最大化原則，使同一節(jié)點任意兩條邊形成的角度盡量大
邊交叉數(shù)量最小原則（最重要）：繪圖時應盡量減少相互交叉邊的數(shù)量。 -
直線邊原則：邊盡量直線，避免曲線，同時線盡可能短。其中，盡量避免交叉是最重要的，然而，并不是所有的圖結構都能滿足這個要求。同時，在性能上，還應當滿足如下需求：
-
算法應當能處理盡可能復雜和龐大的圖，比如上萬甚至上百萬的節(jié)點
- 運算速度還要盡量快，使其能夠在短時間內(nèi)完成操作
- 可交互性，當用戶調(diào)整了某一個點的位置，或改變了一些參數(shù)，算法是否能夠動態(tài)調(diào)節(jié)？
如何表示一個網(wǎng)絡
圖分為有向圖和無向圖，有權圖和無權圖。當圖有權時，權重反映了節(jié)點間的關系，通常值越大，節(jié)點關系越緊密。
- 二維矩陣 最簡單也最方便的表示方法，當網(wǎng)絡結構巨大時，存儲空間浪費是驚人的，而且如果圖是無向圖，則矩陣是對稱的。但存取效率最高
- 臨接表， 即將邊表達為一個數(shù)組，(x,y,dist) *N 對于稀疏圖來說，存儲效率很高，但讀寫效率低。一種簡單的優(yōu)化方法是對x和y 分別進行排序，排序后按二分查找，即可實現(xiàn)快速存取。
- 實時計算 當網(wǎng)絡圖巨大，而計算關系本身的運算并不復雜時，可以采用，此時dist[x][y] 被轉換為一個函數(shù) GetDistance(x,y), 實時完成計算。這是一種很有趣且比較實用的方法，尤其是當網(wǎng)絡圖稀疏時。
為了簡化，我們采用二維矩陣作為存取結構，相信讀者可以很容易的將其改寫為其他形式。
算法介紹
?
力導引算法
力導引布局的方法可以產(chǎn)生相當優(yōu)美的網(wǎng)絡布局，并充分展現(xiàn)網(wǎng)絡的整體結構及其自同構特征。該方法最早由Eades在1984年提出。
基本思想是將網(wǎng)絡看成一個頂點為鋼環(huán)，邊為彈簧的物理系統(tǒng)，系統(tǒng)被賦予某個初始狀態(tài)以后，彈簧彈力（引力和斥力）的作用會導致鋼環(huán)移動，這種運動直到系統(tǒng)總能量減少到最小值停止。每次循環(huán)都要計算任意一對點的斥力和相鄰兩個點的引力，計算復雜度為O（N2）。
基本上絕大多數(shù)算法都遵循著這樣的原則，即：
為了簡化，彈簧的受力并不遵循胡克定律，而是自己定義的受力公式，同時決定引力只在相鄰節(jié)點間。因此，最基本的彈簧算法的偽代碼如下：
隨機分布G的節(jié)點
重復m次
{
計算作用在各個節(jié)點上的力
按照力移動節(jié)點的位置
}
根據(jù)節(jié)點位置繪圖
衡量一個布點算法的美學標準，可以用下面幾個公式來計算：
其中:
?
最后:
改進算法
FR算法改進了彈簧算法，是現(xiàn)在用途最為廣泛的布點算法，很多算法都是在這個算法上改進的。
FR受到了天體重力系統(tǒng)的啟發(fā)，使用力來計算每個節(jié)點的速度，而不是加速度，從而得到每個節(jié)點應當移動的距離。它的每次迭代分為三個步驟：
- 計算節(jié)點之間的排斥力
- 計算相鄰節(jié)點之間的吸引力
- 綜合吸引力和排斥力，通過最大位移限制移動的距離
使用模擬退火算法，使得在圖變得越來越穩(wěn)定時，溫度變得更低，節(jié)點每次移動的距離就變得更小。其主要原因是防止震蕩。
KK算法使得能量最小化，在圖的布局上減少了邊的交叉，除了需要計算所有節(jié)點對之間的最短路徑，并不需要其他理論知識。它雖然每一步的計算復雜度高于FR算法，但迭代次數(shù)較少，使其執(zhí)行速度和效果都比FR好。
實現(xiàn)和性能
上一部分是兩個月前寫的，今天發(fā)現(xiàn)了，想把剩下的補足。
本來想把之前的布點算法的代碼抽出來，做成一個獨立的工程，做成附件發(fā)布出來，奈何之前學弟寫的代碼實在太亂亂亂了，我改了老半天，依然還有一大堆的問題，漸漸地我失去了耐心。覺得這篇文章要太監(jiān)了。
說點感性的認識吧。布點算法很重要，因為能直觀的看出復雜網(wǎng)絡的結構和關系。但在點很多時，性能就得考慮，在時間上，幾百個點的計算幾乎瞬間即可完成，上萬個節(jié)點的圖，采用高性能的布點算法，在犧牲效果的情況下，能在兩分鐘內(nèi)跑完，使用力導引算法，則要花費10分鐘以上。這可能和學弟實現(xiàn)的代碼不夠高效率有關系。
附錄給出了一個簡單的力導引實現(xiàn)，接口應該很容易，一看就能明白。
目前，推薦SM算法，相關鏈接：
https://en.wikipedia.org/wiki/Stress_Majorization
不過，不建議重復造輪子，有現(xiàn)成的graphviz可以使用。應該把時間花在更重要的事情上。如果我當時沒有對語言間互操作那么反感，就不用耽擱那么多功夫重寫C#版本了。
其實原本的代碼中，包含了SM,超大規(guī)模布點，力導引這三種算法的完整代碼，有需要的，可以站內(nèi)我。如果你看到了這里，說明你還是蠻有耐心的。
附錄：C#實現(xiàn)
using System; using System.Collections.Generic; using System.Linq;namespace GraphdotNet.ForceAtlas2Algo {internal class ForceAtlas2{#region Constructors and Destructorspublic ForceAtlas2(int nodeNum){nodes = new Node[nodeNum];for (var i = 0; i < nodes.Length; i++){nodes[i] = new Node();}Degree = new int[nodeNum];}#endregion#region Methodsprivate double GetEdgeWeight(int edgeIdx){return Adjacency[3*edgeIdx + 2];}#endregion#region Constants and Fieldsinternal double Speed;private readonly Node[] nodes;private double outboundAttCompensation = 1;private Region rootRegion;#endregion#region Propertiesinternal bool AdjustSizes { get; set; }internal bool BarnesHutOptimize { get; set; }internal double BarnesHutTheta { get; set; }internal double EdgeWeightInfluence { get; set; }internal double Gravity { get; set; }internal double JitterTolerance { get; set; }internal bool LinLogMode { get; set; }internal bool OutboundAttractionDistribution { get; set; }internal double ScalingRatio { get; set; }internal bool StrongGravityMode { get; set; }private double[] Adjacency { get; set; }private int[] Degree { get; }private int[] EdgeAdjacency { get; set; }#endregion#region Public Methodspublic bool CanAlgo(){return (nodes.Length != 0);}public void EndAlgo(){}public IEnumerable<double> GetX(){return nodes.Select(node => node.X);}public IEnumerable<double> GetY(){return nodes.Select(node => node.Y);}public void GoAlgo(){if (nodes.Length == 0){return;}initialize node datavar nodeNum = nodes.Length;foreach (var n in nodes){n.Olddx = n.Dx;n.Olddy = n.Dy;n.Dx = 0;n.Dy = 0;}// If Barnes Hut active, initialize root RegionImplif (BarnesHutOptimize){rootRegion = new Region(nodes);rootRegion.BuildSubRegions();}// If outboundAttractionDistribution active, compensate.if (OutboundAttractionDistribution){outboundAttCompensation = 0;foreach (var n in nodes){outboundAttCompensation += n.Mass;}outboundAttCompensation /= nodeNum;}// Repulsion (and gravity)var repulsion = ForceFactory.Builder.BuildRepulsion(AdjustSizes, ScalingRatio);const int @from = 0;var to = nodeNum;var rgCalculator = new RepulsionGravityCalculate(nodes,from,to,BarnesHutOptimize,BarnesHutTheta,Gravity,StrongGravityMode ? (ForceFactory.Builder.GetStrongGravity(ScalingRatio)) : (repulsion),ScalingRatio,rootRegion,repulsion);rgCalculator.Run();// Attractionvar attraction = ForceFactory.Builder.BuildAttraction(LinLogMode,OutboundAttractionDistribution,AdjustSizes,1*(OutboundAttractionDistribution ? (outboundAttCompensation) : (1)));var edgeNum = EdgeAdjacency.Length/2;const double epsilon = 1e-6;if (Math.Abs(EdgeWeightInfluence - 0) < epsilon){for (var i = 0; i < edgeNum; i++){attraction.Apply(nodes[EdgeAdjacency[2*i]], nodes[EdgeAdjacency[2*i + 1]], 1);}}else if (Math.Abs(EdgeWeightInfluence - 1) < epsilon){for (var i = 0; i < edgeNum; i++){attraction.Apply(nodes[EdgeAdjacency[2*i]],nodes[EdgeAdjacency[2*i + 1]],GetEdgeWeight(i));}}else{for (var i = 0; i < edgeNum; i++){attraction.Apply(nodes[EdgeAdjacency[2*i]],nodes[EdgeAdjacency[2*i + 1]],Math.Pow(GetEdgeWeight(i), EdgeWeightInfluence));}}// Auto adjust speedvar totalSwinging = 0d; // How much irregular movementvar totalEffectiveTraction = 0d; // Hom much useful movementforeach (var n in nodes){var swinging = Math.Sqrt((n.Olddx - n.Dx)*(n.Olddx - n.Dx) + (n.Olddy - n.Dy)*(n.Olddy - n.Dy));totalSwinging += n.Mass*swinging;// If the node has a burst change of direction, then it's not converging.totalEffectiveTraction += n.Mass*0.5*Math.Sqrt((n.Olddx + n.Dx)*(n.Olddx + n.Dx) + (n.Olddy + n.Dy)*(n.Olddy + n.Dy));}// We want that swingingMovement < tolerance * convergenceMovementvar targetSpeed = JitterTolerance*JitterTolerance*totalEffectiveTraction/totalSwinging;// But the speed shoudn't rise too much too quickly, since it would make the convergence drop dramatically.const double maxRise = 0.5; // Start rise: 50%Speed = Speed + Math.Min(targetSpeed - Speed, maxRise*Speed);// Apply forcesforeach (var n in nodes){// Adaptive auto-speed: the speed of each node is lowered// when the node swings.var swinging = Math.Sqrt((n.Olddx - n.Dx)*(n.Olddx - n.Dx) + (n.Olddy - n.Dy)*(n.Olddy - n.Dy));//double factor = speed / (1f + Math.sqrt(speed * swinging));var factor = Speed/(1 + Speed*Math.Sqrt(swinging));n.X += n.Dx*factor;n.Y += n.Dy*factor;}}public void InitAlgo(double[] adjacencyInput){foreach (var n in nodes){n.ResetNode();}var ran = new Random();foreach (var n in nodes){n.X = (0.01 + ran.NextDouble())*1000 - 500;n.Y = (0.01 + ran.NextDouble())*1000 - 500;}Speed = 1.0;Adjacency = adjacencyInput;var len = Adjacency.Length;var edgeNum = len/3;EdgeAdjacency = new int[edgeNum*2];for (var i = 0; i < edgeNum; i++){EdgeAdjacency[2*i] = (int) Adjacency[3*i];EdgeAdjacency[2*i + 1] = (int) Adjacency[3*i + 1];}for (var i = 0; i < 2*edgeNum; i++){Degree[EdgeAdjacency[i]]++;}var nodeNum = nodes.Length;for (var i = 0; i < nodeNum; i++){nodes[i].Mass = 1 + Degree[i];}}public void OutlierProcess(){var maxRadius = double.NegativeInfinity;foreach (var n in nodes){if (n.Mass >= 2){maxRadius = Math.Max(maxRadius, Math.Sqrt(n.X*n.X + n.Y*n.Y));}}maxRadius += 10;var maxRadiusCeil = (int) Math.Ceiling(maxRadius);var ran = new Random();foreach (var n in nodes){if (n.Mass <= 1){n.X = ran.Next(-maxRadiusCeil*20, maxRadiusCeil*20)/20.0;n.Y = Math.Pow(-1.0, Math.Floor(n.X))*Math.Sqrt(maxRadius*maxRadius - n.X*n.X);}}}public void ResetPropertiesValues(){var nodeNum = nodes.Length;// TuningScalingRatio = nodeNum >= 100 ? 2.0 : 10.0;StrongGravityMode = false;Gravity = 1.0;// BehaviorOutboundAttractionDistribution = false;LinLogMode = false;EdgeWeightInfluence = 1.0;// Performanceif (nodeNum >= 50000){JitterTolerance = 10.0;}else if (nodeNum >= 5000){JitterTolerance = 1.0;}else{JitterTolerance = 0.1;}BarnesHutOptimize = nodeNum >= 1000;BarnesHutTheta = 1.2;}#endregion} }
總結
以上是生活随笔為你收集整理的布点算法的原理和实现的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： Tkinter之Label
- 下一篇： 从URI中获取实际的文件path