我想尽可能有效地模拟网络中随机游走者的运动。下面我展示了一个玩具模型,其中包含我迄今为止尝试过的三种方法。我应该注意到,在我的原始问题中,网络的边缘是固定的,但是边缘的权重可能会更新(即邻居列表是相同的,但权重可能会改变)。
using QuantEcon
using LightGraphs
using Distributions
using StatsBase
n = 700 #number of nodes
#setting an arbitrary network and its transition matrix
G_erdos = erdos_renyi(n, 15/n)
A_erdos = adjacency_matrix(G_erdos) + eye(n, n);
A_transition = A_erdos ./ sum(A_erdos, 2);
##Method 1
#using QuantEcon library
function QE_markov_draw(i::Int, A::Array{Float64,2})
d = DiscreteRV(A[i, :]);
return rand(d, 1)
end
##Method 2
#using a simple random draw
function matrix_draw(i::Int, A::Array{Float64,2}, choices::Array{Int64,1})
return sample(choices, Weights(A[i, :]))
end
##Method 3
# The matrix may be sparse. Therefore I obtain first list of neighbors and weights
#for each node. Then run sample using the list of neighbors and weights.
function neighbor_weight_list(A::Array{Float64,2}, i::Int)
n = size(A)[1]
neighbor_list = Int[]
weight_list = Float64[]
for i = 1:n
for j = 1:n
if A[i, j] > 0
push!(neighbor_list, j)
push!(weight_list, A[i, j])
end
end
end
return neighbor_list, weight_list
end
#Using sample on the reduced list.
function neigh_weights_draw(i::Int, neighs::Array{Int,1}, weigh::Array{Float64,1})
return sample(neighs, Weights(weigh))
end
neighbor_list, weight_list = neighbor_weight_list(A_transition, 1)
states = [i for i = 1:n];
println("Method 1")
@time for t = 1:100000
QE_markov_draw(3, A_transition)
end
println("Method 2")
@time for t = 1:100000
matrix_draw(3, A_transition, states)
end
println("Method 3")
@time for t = 1:100000
neigh_weights_draw(3, neighbor_list, weight_list)
end
一般结果表明(在第一次迭代之后)方法 2 是最快的。方法 3 使用最少的内存,其次是方法 2,但这可能是因为它们在
neighbor_list
上“馈送”和 states
.Method 1
0.327805 seconds (500.00 k allocations: 1.086 GiB, 14.70% gc time)
Method 2
0.227060 seconds (329.47 k allocations: 554.344 MiB, 11.24% gc time)
Method 3
1.224682 seconds (128.19 k allocations: 3.482 MiB)
我想知道哪种实现最有效,以及是否有改进的方法。
最佳答案
以下是我可以给出的一些建议:
在选项 2 中,改用 View 并处理矩阵的转置(因此您处理的是列而不是行):
# here A should be a transpose of your original A
function matrix_draw(i::Int, A::Array{Float64,2}, choices::Array{Int64,1})
return sample(choices, Weights(view(A, i, :)))
end
这在我的测试中提供了近 7 倍的加速。
但一般来说方法 3 应该是最快的,但它似乎没有正确实现。这是一个固定的方法
function neighbor_weight_list(A::Array{Float64,2})
n = size(A)[1]
neighbor_list = Vector{Int}[]
weight_list = Vector{Float64}[]
for i = 1:n
push!(neighbor_list, Int[])
push!(weight_list, Float64[])
for j = 1:n
if A[i, j] > 0
push!(neighbor_list[end], j)
push!(weight_list[end], A[i, j])
end
end
end
return neighbor_list, weight_list
end
function neigh_weights_draw(i::Int, neighs::Vector{Vector{Int}}, weigh::Vector{Vector{Float64}})
return sample(neighs[i], Weights(weigh[i]))
end
neighbor_list, weight_list = neighbor_weight_list(A_transition)
当我运行此代码时,它比固定方法 2 快 4 倍。另请注意,您可以使用方法 3 而不创建邻接矩阵,而是直接使用
neighbors
工作。函数来自 LightGraphs
.
关于julia - 在 Julia 中高效实现马尔可夫链,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/50274934/