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作业君