布点算法的原理和实现
在數據可視化的過程中,繪制網絡拓撲圖是很重要的,它能清晰呈現一個復雜網絡的結構,節點的重要性和關系。比如下面幾張圖:
下面這張圖是我的軟件繪制的:
這些都有一個共同的問題,就是如何讓圖繪制的更加美觀? 復雜的圖結構,已經沒法人工布局了。所以計算機自動布局,就成了一個有趣而重要的話題。我們將其稱為布點算法。
基本原理
一個好的布點算法應該能盡量滿足以下四個特點:
- 對稱性:繪制網絡對應將具有相同結構的子圖圍繞繪圖中心平衡布局
- 連接角最大化原則,使同一節點任意兩條邊形成的角度盡量大
邊交叉數量最小原則(最重要):繪圖時應盡量減少相互交叉邊的數量。 -
直線邊原則:邊盡量直線,避免曲線,同時線盡可能短。其中,盡量避免交叉是最重要的,然而,并不是所有的圖結構都能滿足這個要求。同時,在性能上,還應當滿足如下需求:
-
算法應當能處理盡可能復雜和龐大的圖,比如上萬甚至上百萬的節點
- 運算速度還要盡量快,使其能夠在短時間內完成操作
- 可交互性,當用戶調整了某一個點的位置,或改變了一些參數,算法是否能夠動態調節?
如何表示一個網絡
圖分為有向圖和無向圖,有權圖和無權圖。當圖有權時,權重反映了節點間的關系,通常值越大,節點關系越緊密。
- 二維矩陣 最簡單也最方便的表示方法,當網絡結構巨大時,存儲空間浪費是驚人的,而且如果圖是無向圖,則矩陣是對稱的。但存取效率最高
- 臨接表, 即將邊表達為一個數組,(x,y,dist) *N 對于稀疏圖來說,存儲效率很高,但讀寫效率低。一種簡單的優化方法是對x和y 分別進行排序,排序后按二分查找,即可實現快速存取。
- 實時計算 當網絡圖巨大,而計算關系本身的運算并不復雜時,可以采用,此時dist[x][y] 被轉換為一個函數 GetDistance(x,y), 實時完成計算。這是一種很有趣且比較實用的方法,尤其是當網絡圖稀疏時。
為了簡化,我們采用二維矩陣作為存取結構,相信讀者可以很容易的將其改寫為其他形式。
算法介紹
?
力導引算法
力導引布局的方法可以產生相當優美的網絡布局,并充分展現網絡的整體結構及其自同構特征。該方法最早由Eades在1984年提出。
基本思想是將網絡看成一個頂點為鋼環,邊為彈簧的物理系統,系統被賦予某個初始狀態以后,彈簧彈力(引力和斥力)的作用會導致鋼環移動,這種運動直到系統總能量減少到最小值停止。每次循環都要計算任意一對點的斥力和相鄰兩個點的引力,計算復雜度為O(N2)。
基本上絕大多數算法都遵循著這樣的原則,即:
為了簡化,彈簧的受力并不遵循胡克定律,而是自己定義的受力公式,同時決定引力只在相鄰節點間。因此,最基本的彈簧算法的偽代碼如下:
隨機分布G的節點
重復m次
{
計算作用在各個節點上的力
按照力移動節點的位置
}
根據節點位置繪圖
衡量一個布點算法的美學標準,可以用下面幾個公式來計算:
其中:
?
最后:
改進算法
FR算法改進了彈簧算法,是現在用途最為廣泛的布點算法,很多算法都是在這個算法上改進的。
FR受到了天體重力系統的啟發,使用力來計算每個節點的速度,而不是加速度,從而得到每個節點應當移動的距離。它的每次迭代分為三個步驟:
- 計算節點之間的排斥力
- 計算相鄰節點之間的吸引力
- 綜合吸引力和排斥力,通過最大位移限制移動的距離
使用模擬退火算法,使得在圖變得越來越穩定時,溫度變得更低,節點每次移動的距離就變得更小。其主要原因是防止震蕩。
KK算法使得能量最小化,在圖的布局上減少了邊的交叉,除了需要計算所有節點對之間的最短路徑,并不需要其他理論知識。它雖然每一步的計算復雜度高于FR算法,但迭代次數較少,使其執行速度和效果都比FR好。
實現和性能
上一部分是兩個月前寫的,今天發現了,想把剩下的補足。
本來想把之前的布點算法的代碼抽出來,做成一個獨立的工程,做成附件發布出來,奈何之前學弟寫的代碼實在太亂亂亂了,我改了老半天,依然還有一大堆的問題,漸漸地我失去了耐心。覺得這篇文章要太監了。
說點感性的認識吧。布點算法很重要,因為能直觀的看出復雜網絡的結構和關系。但在點很多時,性能就得考慮,在時間上,幾百個點的計算幾乎瞬間即可完成,上萬個節點的圖,采用高性能的布點算法,在犧牲效果的情況下,能在兩分鐘內跑完,使用力導引算法,則要花費10分鐘以上。這可能和學弟實現的代碼不夠高效率有關系。
附錄給出了一個簡單的力導引實現,接口應該很容易,一看就能明白。
目前,推薦SM算法,相關鏈接:
https://en.wikipedia.org/wiki/Stress_Majorization
不過,不建議重復造輪子,有現成的graphviz可以使用。應該把時間花在更重要的事情上。如果我當時沒有對語言間互操作那么反感,就不用耽擱那么多功夫重寫C#版本了。
其實原本的代碼中,包含了SM,超大規模布點,力導引這三種算法的完整代碼,有需要的,可以站內我。如果你看到了這里,說明你還是蠻有耐心的。
附錄:C#實現
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} }
總結
以上是生活随笔為你收集整理的布点算法的原理和实现的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: Tkinter之Label
- 下一篇: 从URI中获取实际的文件path