Starting the next session five minutes early. So to give our next speaker, Lisa, this
is our own Piscuit policy, no, this Piscuit, hold on, I'm sorry, Piscuit cobble, it's
a tool, okay, apologies. Just to give him five more minutes, because it's a really
depth talk, which is basically velocity sampling efficiently from high dimensional distributions,
I would say you. Okay, thanks a lot. Okay, so today I will try to present this topic,
I will give some brief background and then I will present the software that we have. So let's
start of understanding the title. So what is a distribution, and distribution in let's say
two dimensional spaces, we can have a Gaussian distribution. So for example here the points are
a Gaussian distribution, and we would like to study truncated distribution. So you have points on
the plane that are in Gaussian, and then you have to truncate with some shape. This shape can be
convex or non-convex or linear, non-linear. So here we have just a polygon. So
another distribution, so this is in one dimension, we can have a function, this is a function
that describes the distribution, which is we call it logoncave, and then we want to sample points,
we truncate this distribution, and we want to sample points from minus one to one. And this is
the histograms that we compute. So imagine this objects in with more variables, then we have a
high dimensional distribution, and we want to sample from this. So what is the main, so this is
this are the objects that we want to study, and the problem is that we want to sample from this
distribution. So let's say that you want to compute points, if you have a uniform distribution,
you have a polygon, and you want to compute points. This is a simple example, you want to compute
points that are filling your polygon. So apart from being the fundamental problem in mathematics
and computer science, this sampling is a building block for if you want to do Monte Carlo with
a grayion and volume computation. So why is this interesting, because there are many, many applications,
so for example, in bias and inference, when you want to estimate parameters that have constraints,
so the constraints now is this truncation that we show, your polygon, gives your constraints,
or if you want to do constraint optimization, this is a well known example is a linear programming.
So you have a poly, a poly top, and you want to optimize a function over there.
Another area of application is in finance. In finance, we have, when you have portfolios,
or you play stocks, then you have constraints. So constraints come from some market, or constraints
of constraints over the portfolios, they can be probabilities, and an arbitration is
incapable of biology. So in metabolic networks, when you want to study the metabolism of
organisms, then interestingly, you can model it with the high dimensional poly top. This is
high dimensional example of a polygon, and then the biology is going to sample from there,
from some distribution, and this gives some information about the organism that you want to
study the metabolic. So for finance and biology, I will give certain examples.
I hope that this will be more clear than it's next slide. So if we have a simple,
object, a simple geometric object like a triangle, or,
in fact, the measure we call it a simplex, or we have a hyper cube, then it's easy to sample from
there. We know how to sample. And you can do it also in high dimensional. So they're
formulas for this. You just evaluate the formula and you have your points. But I assume that you have
something, a more strange shape. One way to sample from it is to compute the bonding box,
and then sample from the bonding box, and say, okay, I will sample from the bottom box,
and then see how many points fall inside my object that I want to sample.
Okay, into this is fine, right? So you sample points, all this red and green points,
and say, okay, how many points are the green ones? How many points are all of them? And then
you can estimate the volume of your unknown object, or you can sample from it. It's the same
problem. Okay, this shape, after a few dimensions, so if you draw this in 15 dimensions,
you will not be able to compute the green point. You will sample from the box,
and all points are red. Okay, so in high dimensions we have to do something
and more clever if you want to sample from this unknown shape. Okay, and the solution to this
is random groups. What's a random group? You have a point inside your convex set,
the set that you want to sample from, and then you can move to a neighboring point,
and you continue this procedure for many steps. Okay, so this can be analyzed as a
Markov chain. I will skip all the details about analysis, and if this chain
converted this or not, but just as a general idea, a random group is computing points,
like walks inside your polygon. Okay, and the two examples here is, this is a two
dimensional polygon, but you can do this also in non-convex shapes. Non-linear shapes,
that are convex, you can also do it in non-convex shapes. Okay, so let me now give some examples
of some walks that we use also in practice and these are implemented in our package. So
I will not discuss about the algorithm, I will try to explain it in simple words. So assume that
you learn how to play a billion. When you start playing a billion, you take, you can hit the ball
with not a lot of power, so because you are learning billion. So your ball will be moving inside
the big table, little by little, so you are learning how to play billion. So this is
ten meters of this work, this is the first work, the most simple work, which is called ball work.
So you have a point, you compute a ball with some certain radius, and then you move randomly
in this ball with some probability. If the, if the new point is inside our convex body, we keep
it, otherwise we throw the way. And then you continue likewise in the next step. So you are here,
you compute again a new ball, you sum randomly, and you move. So this is a ball work, this one basic work.
You probably know maybe you know the Metropolis casting algorithm, which is very
has a long history. This is the variant of ball work with our constraints. So if you do this
without the constraints, if you have a point, you compute the ball and you move there, some probability,
this is a Metropolis casting algorithm. Okay, now assume that you are becoming a better player in
billion, and you want to, to suit the balls with more power. Okay, now you are shooting in lines,
so you can suit the ball further away in your table, you are more confident. So this reminds us of
the next work, which is a hidden run. So you are in a point, inside your convex body,
you compute a random line, this randomly, throw randomly from a ball here, a random direction,
this defines a line, you take the segment, which is inside your bolt open, then you move on this segment.
Then you, you draw another line randomly, and you move there on this segment.
So this is a more advanced way to move in your bolt up, because you can go intuitively,
you can imagine that you can go further away. So the work can converge faster, so in practice
converges faster. Okay, and now that we can assume that you are even more advanced, you learn
billion even better, so you want to do something like this. You want to throw your ball,
and then also hit in the borders and do some trajectory reflections. Okay, this is called the
billion work. So you have a point, you select the direction, randomly again, and you suit it there
with some length, so the length is a parameter of this work. So your next point will follow the
trajectories and end up in some other direction. In some other position, inside your
complex body. Okay, this is called the billion work. Okay, just a small comment here that
there are even more advanced works. When where you, you can hit your, you can move inside the
body with some curvature. So here all the, all the trajectories are linear, but you can also
follow some care of the trajectory. And the name of the work is the Hamiltonian model,
the Carlo work, which is based on Hamiltonian dynamics. I will not present it here.
Okay, now one question that comes in mind is how many steps should I do in order to have a point
which is random. So this is an experiment with the works that I present to you. The first row
is the ballwork, the second and the third is the hit and run, and the third is the billion
work, which where you have also reflections. In the beginning we have our work with a very small
step, you can see that the ballwork cannot conferred, but if you give them more steps for a point
to convert, then your sampling is better. So here we try to sum it uniformly from a high
dimensional cube, the cube is 200 dimensions. And this is just a projection in 3B. So hit and run
is a bit better, and interestingly the billion work, even with one step, can convert. So even
one step of a billion work and have a point that is well, you can see it as a random point.
Okay, but in general, let's say that you implement this. It's not possible to draw pictures and see
if our points are nicely distributed. Of course there are theoretical bounds for all these works
or for most of them. That tells you, okay, at the point, if you run this random work for
some bound on the number of the dimensions, the point will conferred, but as you many of you
probably know, in these theoretical bounds, there are huge constants, so they are not practical.
What we do in practice, we can apply some statistical tests, like the effective sum of the size,
or the potential scale reduction factor. So these are test from statistics that you can use to see
if your points are randomly distributed according to the distribution you want to sum.
And a challenge here is that if you use the sum link as a subroutine in some other algorithm that
builds on sampling, then it's not clear how you provide this statistical guarantees.
So these are tools that have theoretical guarantees, those guarantees cannot apply it in practice,
so in practice you do tricks. Okay, so these are these arguments are implemented in a
package called the velocity. It is a C++, and there are two interfaces, one in R and one in Python.
There are different vengeance, VR packets, maybe it's the most mature that it's, it's also,
you can find it in cran, which is the repository for the packages, VR packages, so you can download
it from there and use it. It contains algorithms for sampling, but also, I'll use for volume
integration that builds on sampling, and copul estimation, so copulas are, are these
joint distributions between two or more random variables. So this is how copulas two random
variables look like. This is taking from, from the biologists, this metable with the
biomass of the, of the body. So I'm not explaining it, just to show how copulas look like.
And then it supports optimization for different body, so apart from polytops, which are linear,
you can have a spectahitra, which are, both like this, that appear in semi-defin programming,
and you can have also zonotops, that appear in some applications in robot motion plan.
OK, and the last thing that is contained in this package is utilities for financial and biological
applications, that I will briefly explain later. OK, so now I can briefly show the art in the face
of, of velocity to have an idea how you can use these functions more or less. So here is the,
the example from, something from the square. OK, this is the example of I, I created in this light.
And this is the experiment that I told you that if you go to high dimensions, then you cannot,
you, you cannot compute the semi-correctly. So this is just a, a simple script that computes,
computes points on the square, and then asks if this point falls inside the, the ball,
but goes to, to hire dimensions, as you can see,
into dimensions. So here is the, the dimension, the volume of, the, the, the circle,
the exact volume of the circle and the, the ratio between those two. So you can see that if you
will increase the dimension, up to dimension, this is the, the error starts increasing, and then at some point,
you will not find, in dimension 15, you cannot find the point inside your, uh, certainty. So,
this, rejection sampling, rejection sampling cannot work, if you go to, let's say, more than 10 variables.
OK, so let's see how the, implement the drawdown works are, uh, solving this problem. So,
this is the experiment, uh, with, um, uh, sampling from this high dimensional cube that I,
showing in the slides, uh, here we select, we can select, uh, the work and the work length,
how many steps the algorithm is doing, and, uh, we can have, for example,
uh, if the work step is small, you see that the work is not converges, uh, otherwise for a longer,
uh, step, you have some convergence, and with a long step, uh, you have converges in there for the
bold work, um, this and experiment, we, OK, this is running a sample for, uh, 100 dimensions,
you see that it's quite, uh, fast, uh, so I sample 1000 points in 100 dimensions right now,
uh, with billion work, and this is the result, uh, with only one work step. So,
this converges, uh, really fast. Uh, of course, the cube is not, it's, it's a simple, uh,
convex region, but, uh, it's good for, uh, the demonstration.
OK, and then this is the, uh, I can run all of them, and, uh, uh, compute some statistical tests,
so this is the ESS test. It tells you more or less how many points are, um,
how many points are statistical independent. So, let's say here, we compute 1000 points in 100
dimensions, uh, with hidden run, uh, only three or 19 are independent. Let's say that random
is distributed nicely. Uh, in bulk work is even worse, and in billiard work, as we also saw,
uh, before you have a very big number. So more than half of the points are,
nicely distributed. So this gives you a quantitative, uh, test of the pictures that we show before.
So we can use it in practice. So we can implement, uh, an algorithm and, uh, computer
that's a, this is called test, and if this statistical test is, uh, good, we, we continue. So we can
actually use it. It's not just like a plot. Uh, OK, and then, uh, we can, I, I'm not describing the
algorithm, but using sample, you can, uh, you can implement, uh, algorithms for, uh, volume computation.
So let's say that we want to compute the volume of, um, of a cube in four dimensions.
We know that this will be 16. Uh, so this randomize algorithms. If I land many times, I get different
values around around 16. So I can do, uh, statistics there that can increase the, um,
varied bound. Um, but I only want to show this. Uh, so here I can, I can, I can, I can compute statistics
on the, I can compute many volumes and then do some, uh, mean and median, uh, uh,
computations to have some statistics. Uh, I only want to show one example here that shows the,
the power of these methods. Um, here is a well known, uh, a portal from that,
combinatorics. This is very difficult to compute the volume. Um, the 10th version of this
polytope is, uh, in 81 dimensions and the volume is this number, this fraction. This is computed
by a massive parallel computation after, um, I think weeks. Um, and, uh, this is the exact volume.
So we can run our software on this polytope to have an estimation. And this, uh, only takes,
some, uh, seconds. And we have something that hits the order of magnitude. Of course, we
kind of expect to have, uh, the same number. Uh, but, uh, we can, uh, we can use it for, uh, up to 1000
dimensions and compute the volume of these polytopes in 1000 dimensions. So, this, so, uh,
showcase the, the scaling of, uh, these methods. And, uh, you can also use it in the
regression. So you have a polynomial and you want to integrate it over, uh, a convict set.
Uh, this can be either, this can be done by this packet screw back to, but since this is deterministic,
it cannot go to high dimensions. Uh, I have an example here, uh, with the list that can compute,
uh, a good estimation, uh, a good estimation of, uh, this is Degral. And there is an interesting
Python, uh, packets by, uh, European space agency that computes, uh, this kind of integrals,
uh, with the GPUs. So I didn't test, uh, the comparison, but, uh, they, they try to, to, to overcome
the, the problem of, uh, uh, uh, uh, high dimensions with the GPUs. So this also,
any interesting, uh, direction. Okay, let me go back to slides to slides.
And just so the two applications that, uh, uh, we can have, uh, with sampling. So first, uh,
finance allocation. So the set of portfolios, uh, with this investment on a set of stocks,
uh, can be represented as a simplex. So, uh, if you have two stocks, then, uh, you have a triangle.
And all the possible portfolios are points in this triangle. Okay. So you can also have extra
constraints, uh, so you take this triangle and you cut it by new constraints like, okay, I want to play
stocks, only, uh, U.S., or whatever. And you create some polygon. So imagine this, if you
don't have only two stocks, if you have thousands of stocks, then this is a huge, uh, uh, political.
And you want to study this, uh, we started this particular scenario where portfolios have the
same volatility. So volatility is the degrees of a variation of, uh, of the prize. So
this, this means that all the stocks with the same volatility, they are, are, they lay on the surface of
a ellipseoid. So into these, uh, is an ellipse. Uh, so you want an intersection of these ellipse
with a triangle. Uh, as you can imagine in, even in three dimensions, this shape is crazy. It's non-convex.
Uh, and, uh, we, we implemented some, uh, algorithms based of these sampling
angles that I explained you before with random walks, uh, that started, uh, uh, the anomaly detection,
in those stock markets. So the talent is to, to sum up from this, uh, objects that are not
complex and in high dimensions. Okay. And, uh, an application in structure biology,
interestingly, again, uh, if you go to, to study with a volcanic networks, uh, again, you have, uh,
a polytop or into the, you have a polygon. And the points, uh, there are, uh, um,
specific cases of, uh, uh, of the flow in your network. Let's say,
kind of adding that there is a flow in the metabolic network, because, uh, they, they are,
connections between, uh, metabolites. Uh, so, uh, you, you want, you, if, if the biologist
is going to study the balance state, uh, the state state of, uh, this, um, uh, this system.
And this is a commercial job. So you can use exactly the same, uh, algorithms there, uh,
and, uh, you, you can study, uh, for example, uh, uh, uh, uh, uh, you, you can start the
some reactions and, uh, compute histograms like this, uh, by sampling. So, this is an interesting
challenge because those polytop are, uh, since the networks, those metabolic networks can be
very large, uh, are in thousands of dimensions. So, there are not many alternatives, uh,
up to now to, to study those networks with, uh, sampling. So you have to use some, some,
implementation like this. Okay, uh, just to wrap up, uh, this, uh, these algorithms for sampling
and volume, uh, integration and implemented the C++. It's, uh, they're open source. So,
germ scale is just the set, uh, in GitHub of, uh, all these, uh, packages, uh, it's the C++
package and the two interfaces. Uh, yes, I like these, uh, QR that the previous talks, uh,
had, so I also had the QR. So, yes, this is supported by, uh, open source community
and people that volunteer to, if, if they have, they find the interesting to, to work on this.
Uh, so, this basic, basically how this is implemented and maintained. Okay, so, that's all.
Thank you.