辅导案例-PHY 451/551

欢迎使用51辅导,51作业君孵化低价透明的学长辅导平台,服务保持优质,平均费用压低50%以上! 51fudao.top
Bayesian Data Analysis
A PHY 451/551 and I CSI 451/551
HW#9p

Your Mission this week is to:
Code Skilling’s Nested Sampling Algorithm into Matlab


This is facilitated by the fact that the code is already written.
Although it is made difficult by the fact that it is written in C and not Matlab.

Here are a few hints to help guide you through the process:

This will be a little difficult at first, but it will help you to better learn to read programs
and modify them. Be sure to take the time to read through Skilling’s code and try to
imagine what it is doing. In addition read through the code I have attached, compare it to
Skilling’s code to see how it is similar and different.

Test parts of the code separately to make sure it works properly.

Pay attention to the messages that the errors give you…as annoying as they are, they are
trying to help.

Read chapter 9.
It is not an easy read, but do the best you can.

Section 9.2.4 introduces the main nested sampling code. This code is written as a
C function called main.c. We will be calling it main.m. The idea is that this
function must call several other functions that are application specific. The
benefit is that main.m will be able to be used on a wide variety of problems with
little-to-no modification.

Section 9.3.2 introduces the rest of the code in terms of an application called The
Lighthouse Problem. You read about this problem in Chapter 2 although we have
not yet discussed it. There are 5 functions in this section:
apply
logLhood
Prior
Explore
Results
which are all included under apply.c along with a structure definition.

The code is written as functions…not a script!
Write:
main.m
apply.m
logLhood.m
Prior.m
Explore.m
Results.m
as separate m-file functions.

Start with the easy one first…and go from there.


There are significant differences between C and Matlab

First is the typedef struct { } object; structure
In C this is used to make a structure which stores variables conveniently.
Think of it as a box.
We will use a Cell Array, which you have seen before.

I set this up in apply.m, which I have included to help you get started.
Read it over closely and note the similarities and differences. I define the Cell Arrays in
apply.m, and I define the data here as well (instead of in logLhood.m like Skilling does).

Once the cell arrays are defined, you can access their contents easily.
I have attached prior.m as well, which uses the cell arrays in Matlab rather than the C
structures. Read my code and Skilling’s code and note the similarities and differences.

Also in Prior.m, I use the function rand() or just rand, to generate a random uniformly
distributed number. This is in contrast to the #define statement that Skilling uses in
main.c to define UNIFORM as a random number generated between 0 and 1.

Skilling also defines a logarithmic addition in main.c called PLUS(x,y).
He uses it later on in a line that says:
logZnew = PLUS(logz, Obj[worst].logWt);
You should implement it this way:

% implement PLUS here
if (logZ > Obj(worst).logWt)
logZnew = logZ + log(1 + exp(Obj(worst).logWt-logZ));
else
logZnew = Obj(worst).logWt + log(1 + exp(logZ –Obj(worst).logWt));
end

Last, there are a few other bumps:
- Comments in C are //
- Be aware that in C arrays are accessed by Obj[worst] whereas in Matlab we use
Obj(worst).
- When he calls Prior and Explore from main, he does so with the addresses of the
objects &Obj… You can just pass the object cell array, as in:

% Set prior objects
for i = 1:n
Obj(i) = Prior(Obj(i));
end

% Evolve copied object within constraint
Obj(worst) = Explore(Obj(worst), Try, logLstar);

% Kill Worst Object in favor of a copy of a different survivor
copy = ceil(n * rand()); % choose an object between 1 and n
while ((copy == worst) && n>1)
copy = ceil(n * rand()); % choose another object
end
logLstar = Obj(worst).logL; % new likelihood constraint
Obj(worst) = Obj(copy); % overwrite worst object


This will be a little difficult at first, but it will help you to better learn to read
programs and modify them.


欢迎咨询51作业君
51作业君

Email:51zuoyejun

@gmail.com

添加客服微信: abby12468