@@ -11,8 +11,12 @@ include("lib_pbp.jl")
11
11
include (" lib_ep.jl" )
12
12
include (" lib_doSim.jl" )
13
13
#
14
- demoNames = [" demoGrid" ," demoChain" ," demoImg" ]
15
- expname = demoNames[1 ] # choice of demo (reason for this syntax is to show possibilities)
14
+ demoNames = [" demoGrid" ," demoTree" ," demoImg" ," demoChain" ]
15
+ expname = demoNames[2 ] # choice of demo (reason for this syntax is to show possibilities)
16
+ # make directory to store stuff
17
+ if ~ isdir (expname)
18
+ mkdir (expname)
19
+ end
16
20
#
17
21
# ==================================================================================================
18
22
if expname == " demoGrid"
@@ -83,77 +87,72 @@ if expname == "demoGrid"
83
87
s_init = 4 * obs_var
84
88
end
85
89
# ==================================================================================================
86
- elseif expname == " demoChain "
90
+ elseif expname == " demoTree "
87
91
#
88
92
# SIMULATIONS TO BE RUN
89
93
#
90
94
RELOAD = true # re-generate everything
91
95
LBPD = true # LBP on deterministic grid
92
96
EPBP = true # EPBP
93
97
FEPBP = false # Fast-EPBP
94
- PBP = false # PBP with MH sampling
95
- EP = true # straight EP
98
+ PBP = true # PBP with MH sampling
99
+ EP = false # straight EP
96
100
#
97
101
# VERBOSITY
98
102
NODE_PROGRESS = false
99
103
#
100
104
# SIMULATION PARAMETERS [!USER!]
101
105
#
102
- Nlist = [100 , 200 ] # (list) number of particles per node
103
- Clist = [7 , 10 ] # (list) number of components for FEPBP, need to be of same dim as NLIST
106
+ Nlist = [100 ] # (list) number of particles per node
107
+ Clist = [7 ] # (list) number of components for FEPBP, need to be of same dim as NLIST
104
108
Ninteg = 30 # number of integration points for EP proj
105
109
Ngrid = 200 # number of points in the discretization
106
110
nloops = 10 # number of loops through scheduling
107
111
nEPloops = 25 # number of EP iterations
108
112
nruns = 1 # number of time we run the whole thing
109
113
#
110
- # EP PROJECTION MODE, default is KL ignoring update if moments not valid.
114
+ # EP PROJECTION MODE, default is KL ignoring update if moments not valid
111
115
EP_PROJ_MLE = false # use MLE projection instead of KL-EP
112
116
#
113
117
# Additional parameters for PBP
114
118
#
115
119
MHIter = 20 # number of MH iterations
116
- MHProposal = Normal (0 ,. 1 ) # form of the MH proposal
120
+ MHProposal = Normal (0 ,1 ) # form of the MH proposal
117
121
PARACHAINS = true
118
122
#
119
123
# DECLARE GM
120
124
#
121
- T = 5
122
- nnodes,nedges,edge_list = gm_chain (T)
123
- # > declare scheduling
124
- scheduling = gm_chain_scheduling (T,true ) # forward only
125
+ nnodes,nedges = 8 ,7
126
+ edge_list = [ 1 2 ; 1 3 ; # top 1-{2,3}
127
+ 2 1 ; 2 4 ; 2 5 ; 2 6 ; 3 1 ; 3 7 ; 3 8 ; # middle 2-{1,4,5,6} and 3-{1,7,8}
128
+ 4 2 ; 5 2 ; 6 2 ; 7 3 ; 8 3 ; ] # leaves 4-2,5-2,6-2,7-3,8-3
129
+ sched1 = hcat (1 ,2 ,4 ,5 ,6 ,3 ,7 ,8 )
130
+ sched2 = hcat (4 ,5 ,6 ,2 ,7 ,8 ,3 ,1 )
131
+ scheduling = vcat (sched1[:],sched2[:])
125
132
# > declare edge and node potential
126
- HOMOG_EDGE_POT = false # if edge pot is symmetric & the same everywhere (eg: image)
127
- #
128
- node_noise = Normal (0 ,2 )
129
- node_mult = 0.7
130
- edge_noise = Normal (0 ,1 )
131
- edge_mult = 0.5
132
- #
133
- eval_edge_pot (from,to,xfrom,xto) = pdf (edge_noise,xto- edge_mult* xfrom)
134
- eval_node_pot (node,xnode) = pdf (node_noise,obs_values[node]- node_mult* xnode)
135
- #
136
- # > sampling from MH?
133
+ HOMOG_EDGE_POT = true # edge pot is symmetric & same everywhere (eg: image)
134
+ # > side functions to define edge/node potential
135
+ node_potential = MixtureModel ([Normal (- 2 ,1 ),Gumbel (1 ,1.5 )],[0.7 ,0.3 ])
136
+ edge_potential = Cauchy (0 ,1 )
137
+ # > definition of edge/node potential functions
138
+ eval_edge_pot (from,to,xfrom,xto) = pdf (edge_potential,xfrom- xto)
139
+ eval_node_pot (node,xnode) = pdf (node_potential,obs_values[node]- xnode)
140
+ # > (PBP) sampling from MH?
137
141
sampleMHP (old) = old+ rand (MHProposal,N)'
142
+ LMHCHAIN = 1000
143
+ ENDCHUNK = 300
144
+ sampleMHP2 (old) = old+ rand (Normal (0 ,0.1 ),1 )
145
+ # > initial values on the graph
146
+ orig_values = zeros (nnodes,1 ) + 2
138
147
#
139
- est_range = (- 10 , 10 ) # > estimated 1D-range for integration
148
+ est_range = (- 7 , 7 ) # > estimated 1D-range for integration
140
149
sigma_thresh = 0.01
141
- # > generate observations
142
150
if RELOAD
143
- orig_values = zeros (nnodes,1 )
144
- obs_values = zeros (nnodes,1 )
145
- #
146
- orig_values[1 ] = rand (Normal (0 ,1 ))
147
- obs_values[1 ] = rand (Normal (node_mult* orig_values[1 ],1 ))
148
- #
149
- for node= 2 : T
150
- orig_values[node] = rand (Normal (edge_mult* orig_values[node- 1 ],2 ))
151
- obs_values[node] = rand (Normal (node_mult* orig_values[node],1 ))
152
- end
151
+ # > generate observations
152
+ obs_values = orig_values + rand (node_potential,nnodes)
153
153
writecsv (" $expname /$expname \_ orig_values.dat" ,orig_values)
154
154
writecsv (" $expname /$expname \_ obs_values.dat" ,obs_values)
155
155
#
156
- # > to start
157
156
obs_var = sqrt (var (obs_values))
158
157
s_init = 4 * obs_var
159
158
end
@@ -165,9 +164,9 @@ elseif expname == "demoImg"
165
164
RELOAD = false # re-generate everything
166
165
LBPD = false # LBP on deterministic grid
167
166
EPBP = false # EPBP
168
- FEPBP = true # Fast-EPBP
167
+ FEPBP = true # Fast-EPBP
169
168
PBP = false # PBP with MH sampling
170
- EP = false # straight EP
169
+ EP = false # straight EP
171
170
#
172
171
# VERBOSITY
173
172
NODE_PROGRESS = false
@@ -236,16 +235,89 @@ elseif expname == "demoImg"
236
235
function eval_node_pot (node,xnode)
237
236
return pdf (node_potential,obs_values[node]- xnode)
238
237
end
238
+ # ==================================================================================================
239
+ elseif expname == " demoChain"
240
+ #
241
+ # SIMULATIONS TO BE RUN
242
+ #
243
+ RELOAD = true # re-generate everything
244
+ LBPD = true # LBP on deterministic grid
245
+ EPBP = true # EPBP
246
+ FEPBP = false # Fast-EPBP
247
+ PBP = false # PBP with MH sampling
248
+ EP = true # straight EP
249
+ #
250
+ # VERBOSITY
251
+ NODE_PROGRESS = false
252
+ #
253
+ # SIMULATION PARAMETERS [!USER!]
254
+ #
255
+ Nlist = [100 ,200 ] # (list) number of particles per node
256
+ Clist = [7 ,10 ] # (list) number of components for FEPBP, need to be of same dim as NLIST
257
+ Ninteg = 30 # number of integration points for EP proj
258
+ Ngrid = 200 # number of points in the discretization
259
+ nloops = 10 # number of loops through scheduling
260
+ nEPloops = 25 # number of EP iterations
261
+ nruns = 1 # number of time we run the whole thing
262
+ #
263
+ # EP PROJECTION MODE, default is KL ignoring update if moments not valid.
264
+ EP_PROJ_MLE = false # use MLE projection instead of KL-EP
265
+ #
266
+ # Additional parameters for PBP
267
+ #
268
+ MHIter = 20 # number of MH iterations
269
+ MHProposal = Normal (0 ,.1 ) # form of the MH proposal
270
+ PARACHAINS = true
271
+ #
272
+ # DECLARE GM
273
+ #
274
+ T = 5
275
+ nnodes,nedges,edge_list = gm_chain (T)
276
+ # > declare scheduling
277
+ scheduling = gm_chain_scheduling (T,true ) # forward only
278
+ # > declare edge and node potential
279
+ HOMOG_EDGE_POT = false # if edge pot is symmetric & the same everywhere (eg: image)
280
+ #
281
+ node_noise = Normal (0 ,2 )
282
+ node_mult = 0.7
283
+ edge_noise = Normal (0 ,1 )
284
+ edge_mult = 0.5
285
+ #
286
+ eval_edge_pot (from,to,xfrom,xto) = pdf (edge_noise,xto- edge_mult* xfrom)
287
+ eval_node_pot (node,xnode) = pdf (node_noise,obs_values[node]- node_mult* xnode)
288
+ #
289
+ # > (PBP) sampling from MH?
290
+ sampleMHP (old) = old+ rand (MHProposal,N)'
291
+ LMHCHAIN = 1000
292
+ ENDCHUNK = 300
293
+ sampleMHP2 (old) = old+ rand (Normal (0 ,0.1 ),1 )
294
+ #
295
+ est_range = (- 10 ,10 ) # > estimated 1D-range for integration
296
+ sigma_thresh = 0.01
297
+ # > generate observations
298
+ if RELOAD
299
+ orig_values = zeros (nnodes,1 )
300
+ obs_values = zeros (nnodes,1 )
301
+ #
302
+ orig_values[1 ] = rand (Normal (0 ,1 ))
303
+ obs_values[1 ] = rand (Normal (node_mult* orig_values[1 ],1 ))
304
+ #
305
+ for node= 2 : T
306
+ orig_values[node] = rand (Normal (edge_mult* orig_values[node- 1 ],2 ))
307
+ obs_values[node] = rand (Normal (node_mult* orig_values[node],1 ))
308
+ end
309
+ writecsv (" $expname /$expname \_ orig_values.dat" ,orig_values)
310
+ writecsv (" $expname /$expname \_ obs_values.dat" ,obs_values)
311
+ #
312
+ # > to start
313
+ obs_var = sqrt (var (obs_values))
314
+ s_init = 4 * obs_var
315
+ end
239
316
end
240
317
241
318
#
242
319
# ======== RUN SIMULATIONS =========================================================================
243
320
#
244
- # make directory to store stuff
245
- if ~ isdir (expname)
246
- mkdir (expname)
247
- end
248
- #
249
321
integ_pts = linspace (est_range[1 ],est_range[2 ],Ninteg)' # ! leave the transpose
250
322
grid = linspace (est_range[1 ],est_range[2 ],Ngrid)'
251
323
# > generate observations
0 commit comments