Arquivos da categoria: Hardcore

Para matemáticos e físicos que passarem por aqui.

Propagating randomness

Hardcore

Ricardo: O post de hoje é em inglês, e de um convidado muito especial. Alexander Dobrinevski, autor do blog inordinatum.wordpress.com, é meu grande amigo bielo-russo que, após sua graduação na Universidade Ludwig-Maximilians de Munique, seguiu para dois anos de mestrado em estudos avançados em matemática na Universidade de Cambridge e deu-me o privilégio de sua companhia no Laboratório de Física Teórica da École Normale Supérieure, onde dividimos sala durante alguns meses de trabalho enquanto ele fazia seu doutorado e eu, o mestrado. Viciado em cafés, filmes não tão mainstream e transformadas de Laplace. O post de hoje é bem hardcore, vale aos especialistas da área, e a mim, que aprendi e aprendo bastante a cada vez que nos falamos.

Ricardo: This post will be in english, for we have a very special guest today. Alexander Dobrinevski, author of the blog inordinatum.wordpress.com, is a good friend of mine from Belarus. Having finished his graduation at LMU, he followed the Advanced Studies in Mathematics program from Cambridge University to end up his doctoral thesis at the École Normale, where I had the pleasure to share with him the office for a few months. He is addicted to coffee, exotic movies and Laplace transforms. Today’s post is also very hardcore, to the benefit of our dear specialist in the field, and to mine, I who have learned so much and still do every time I get to talk to this dearest friend of mine.


Introduction

I guess anybody doing statistical physics or probability theory has played around with Brownian Motion (by which I mean, here and in the following, the Wiener process, and not the physical phenomenon) at some time or another. It is used e.g. for modelling the price of a stock, the position of a particle diffusing in a gas or liquid, or the pinning force on an elastic interface in a disordered medium.

Being Markovian, time evolution of Brownian motion is completely determined by its propagator, i.e. the probability (density) to arrive at $x_1$ at time $ t_1$ starting from $ x_0$ at time $ t_0$. This is, of course, known to be a Gaussian. However, for practical applications, one often needs to restrict the Wiener process (in general, with drift) to a half-line or an interval, with some imposed boundary conditions (absorbing or reflecting). In the example of a stock price, these would describe call or put options on the stock. In this post I will derive, hopefully in a pedagogical way, the propagator of Brownian motion with drift and linear boundaries.

Propagators without drift

Let us first consider standard Brownian motion $ W(t)$ without drift. Without loss of generality, let us assume it starts at $W(0)=0$. It satisfies the Langevin equation

\[\dot{W}(t) = \xi(t)\]

where $ \xi(t)$ is Gaussian white noise with correlation

\[\overline{\xi(t)\xi(t’)} = 2\sigma \delta(t-t’).\]

With these conventions, the free propagator (i.e. the propagator without any boundaries) $ P(x,t)$ is given by the solution of the Fokker-Planck equation

\[ \partial_t P(x,t) = \sigma \partial_x^2 P(x,t)\]

with initial condition $ P(x,0) = \delta(x)$. This PDE, also known as the heat equation, is easily solved by taking a Fourier transform. The solution is given by

\[ P(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}e^{-\frac{x^2}{4\sigma t}}.\]

Now, let us determine the propagator with an absorbing boundary at $ x=b > 0$. In the Fokker-Planck equation, this is equivalent to the boundary condition $ P(b,t)=0$, which makes applying Fourier transforms difficult. However, we can use the method of images to find the solution: $ P(b,t)=0$ is enforced automatically if we add a negative source at $ x=2b$ (the position of the original source, reflected at $ b$), i.e. take the initial condition $ P(x,0) = \delta(x)-\delta(x-2b)$. The final propagator with an absorbing boundary at $ x=b$ is thus

\[P^{(b)}(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}\left[e^{-\frac{x^2}{4\sigma t}}-e^{-\frac{(x-2b)^2}{4\sigma t}}\right].\]

Similarly one can treat the case of two absorbing boundaries, one at $ x=b>0$, one at $ x=a<0$. One then needs an infinite series of images, and obtains the propagator as a series which can be rewritten in terms of Jacobi Theta functions.

Propagators with drift

Now let us generalize to the Brownian motion with drift $ \mu$. Then the Langevin equation for $ W(t)$ becomes

\[ \dot{W}(t) =\ mu + \xi(t).\]

The free propagator is obtained from the Fokker-Planck equation just as above:

\[P(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}e^{-\frac{(x-\mu t)^2}{4\sigma t}}.\]

Let us now introduce again a constant absorbing boundary at $ x=b$. Applying the method of images is not so straightforward anymore. Due to the drift, a path which goes from $ x=0$ to the boundary $ x=b$ will not have the same weight as the reflected path which goes from $ x=2b$ to the boundary $ x=b$. However, for the case of constant drift considered here, the weights of Brownian paths with and without drift have a simple relationship. In my view, the easiest way to see it is using path integrals. The propagator is given by

\[ P(x_f,t) = \int_{x(0)=0}^{x(t)=x_f} \mathcal{D}[x]e^{-\int_0^t \mathrm{d}s, \frac{1}{4\sigma}\left(\dot{x}(s)-\mu\right)^2}\]

Now, expanding the “action” in the exponent, using the fact that our drift $ \mu$ is constant, and using our boundary conditions, this is equal to

\[ P(x_f,t) = e^{-\frac{\mu}{2\sigma}x_f+ \frac{\mu^2}{4\sigma}t}\int_{x(0)=0}^{x(t)=x} \mathcal{D}[x]e^{-\int_0^t \mathrm{d}s, \frac{1}{2\sigma}\left(\dot{x}(s)\right)^2}.\]

We thus get a simple weight depending on the final position, but the remaining path integral is taken over a drift-less Brownian motion, and there we know the solution already, both with and without the boundary! In mathematical literature, you will often find this manipulation under the name of the Cameron-Martin-Girsanov theorem, but I find the path integral explanation much clearer for somebody coming from physics. Note that in the case where the drift $ \mu$ is a function of time, we cannot pull the weight out of the path integral, because it involves the whole trajectory and not just the final point. This shows why non-constant drift with absorbing boundaries is a much more complicated problem (although the free propagator is still trivial to write down!).

The final formula for the propagator of the Brownian motion with drift $ \mu$ and an absorbing boundary at $ x=b$ (also known as Bachelier-Levy formula) is thus

\[P^{(b,mu)}(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}e^{-\frac{mu}{2\sigma}x+ \frac{mu^2}{4\sigma}t}\left[e^{-\frac{x^2}{4\sigma t}}-e^{-\frac{(x-2b)^2}{4\sigma t}}\right]\]

This is now all that is required to compute things like first-passage times, survival probabilities, etc. The generalization to two absorbing boundaries follows from the solution in the driftless case by multiplying with the same weight as here.

I hope you see that with the right tools, obtaining these propagators is nothing miraculous. If you have any questions or comments, I’d be glad to hear them! If people are interested, at some point I may write a continuation of this blog post, possibly on generalizations of the methods discussed here to Ornstein-Uhlenbeck processes, Bessel Processes, or more complicated boundaries.

Have fun, and thanks for reading!

Ricardo: This post will be in english, for we have a very special guest today. Alexander Dobrinevski, author of the blog inordinatum.wordpress.com, is a good friend of mine from Belarus. Having finished his graduation at LMU, he followed the Advanced Studies in Mathematics program from Cambridge University to end up his doctoral thesis at the École Normale, where I had the pleasure to share with him the office for a few months. He is addicted to coffee, exotic movies and Laplace transforms. Today’s post is also very hardcore, to the benefit of our dear specialist in the field, and to mine, I who have learned so much and still do every time I get to talk to this dearest friend of mine.


Introduction

I guess anybody doing statistical physics or probability theory has played around with Brownian Motion (by which I mean, here and in the following, the Wiener process, and not the physical phenomenon) at some time or another. It is used e.g. for modelling the price of a stock, the position of a particle diffusing in a gas or liquid, or the pinning force on an elastic interface in a disordered medium.

Being Markovian, time evolution of Brownian motion is completely determined by its propagator, i.e. the probability (density) to arrive at $x_1$ at time $ t_1$ starting from $ x_0$ at time $ t_0$. This is, of course, known to be a Gaussian. However, for practical applications, one often needs to restrict the Wiener process (in general, with drift) to a half-line or an interval, with some imposed boundary conditions (absorbing or reflecting). In the example of a stock price, these would describe call or put options on the stock. In this post I will derive, hopefully in a pedagogical way, the propagator of Brownian motion with drift and linear boundaries.

Propagators without drift

Let us first consider standard Brownian motion $ W(t)$ without drift. Without loss of generality, let us assume it starts at $W(0)=0$. It satisfies the Langevin equation

\[\dot{W}(t) = \xi(t)\]

where $ \xi(t)$ is Gaussian white noise with correlation

\[\overline{\xi(t)\xi(t’)} = 2\sigma \delta(t-t’).\]

With these conventions, the free propagator (i.e. the propagator without any boundaries) $ P(x,t)$ is given by the solution of the Fokker-Planck equation

\[ \partial_t P(x,t) = \sigma \partial_x^2 P(x,t)\]

with initial condition $ P(x,0) = \delta(x)$. This PDE, also known as the heat equation, is easily solved by taking a Fourier transform. The solution is given by

\[ P(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}e^{-\frac{x^2}{4\sigma t}}.\]

Now, let us determine the propagator with an absorbing boundary at $ x=b > 0$. In the Fokker-Planck equation, this is equivalent to the boundary condition $ P(b,t)=0$, which makes applying Fourier transforms difficult. However, we can use the method of images to find the solution: $ P(b,t)=0$ is enforced automatically if we add a negative source at $ x=2b$ (the position of the original source, reflected at $ b$), i.e. take the initial condition $ P(x,0) = \delta(x)-\delta(x-2b)$. The final propagator with an absorbing boundary at $ x=b$ is thus

\[P^{(b)}(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}\left[e^{-\frac{x^2}{4\sigma t}}-e^{-\frac{(x-2b)^2}{4\sigma t}}\right].\]

Similarly one can treat the case of two absorbing boundaries, one at $ x=b>0$, one at $ x=a<0$. One then needs an infinite series of images, and obtains the propagator as a series which can be rewritten in terms of Jacobi Theta functions.

Propagators with drift

Now let us generalize to the Brownian motion with drift $ \mu$. Then the Langevin equation for $ W(t)$ becomes

\[ \dot{W}(t) =\ mu + \xi(t).\]

The free propagator is obtained from the Fokker-Planck equation just as above:

\[P(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}e^{-\frac{(x-\mu t)^2}{4\sigma t}}.\]

Let us now introduce again a constant absorbing boundary at $ x=b$. Applying the method of images is not so straightforward anymore. Due to the drift, a path which goes from $ x=0$ to the boundary $ x=b$ will not have the same weight as the reflected path which goes from $ x=2b$ to the boundary $ x=b$. However, for the case of constant drift considered here, the weights of Brownian paths with and without drift have a simple relationship. In my view, the easiest way to see it is using path integrals. The propagator is given by

\[ P(x_f,t) = \int_{x(0)=0}^{x(t)=x_f} \mathcal{D}[x]e^{-\int_0^t \mathrm{d}s, \frac{1}{4\sigma}\left(\dot{x}(s)-\mu\right)^2}\]

Now, expanding the “action” in the exponent, using the fact that our drift $ \mu$ is constant, and using our boundary conditions, this is equal to

\[ P(x_f,t) = e^{-\frac{\mu}{2\sigma}x_f+ \frac{\mu^2}{4\sigma}t}\int_{x(0)=0}^{x(t)=x} \mathcal{D}[x]e^{-\int_0^t \mathrm{d}s, \frac{1}{2\sigma}\left(\dot{x}(s)\right)^2}.\]

We thus get a simple weight depending on the final position, but the remaining path integral is taken over a drift-less Brownian motion, and there we know the solution already, both with and without the boundary! In mathematical literature, you will often find this manipulation under the name of the Cameron-Martin-Girsanov theorem, but I find the path integral explanation much clearer for somebody coming from physics. Note that in the case where the drift $ \mu$ is a function of time, we cannot pull the weight out of the path integral, because it involves the whole trajectory and not just the final point. This shows why non-constant drift with absorbing boundaries is a much more complicated problem (although the free propagator is still trivial to write down!).

The final formula for the propagator of the Brownian motion with drift $ \mu$ and an absorbing boundary at $ x=b$ (also known as Bachelier-Levy formula) is thus

\[P^{(b,mu)}(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}e^{-\frac{mu}{2\sigma}x+ \frac{mu^2}{4\sigma}t}\left[e^{-\frac{x^2}{4\sigma t}}-e^{-\frac{(x-2b)^2}{4\sigma t}}\right]\]

This is now all that is required to compute things like first-passage times, survival probabilities, etc. The generalization to two absorbing boundaries follows from the solution in the driftless case by multiplying with the same weight as here.

I hope you see that with the right tools, obtaining these propagators is nothing miraculous. If you have any questions or comments, I’d be glad to hear them! If people are interested, at some point I may write a continuation of this blog post, possibly on generalizations of the methods discussed here to Ornstein-Uhlenbeck processes, Bessel Processes, or more complicated boundaries.

Have fun, and thanks for reading!

Introduction

Quiconque ayant déjà pratiqué la physique statistique ou les probabilités a forcément eu affaire au mouvement Brownien (j’entends par là le processus de Wiener, et pas le phénomène physique). On l’utilise, par exemple, pour modéliser l’évolution des prix, la position des particules dans un liquide ou un gaz, ou la force d’ancrage sur l’interface élastique d’un milieu désordonné.

L’évolution temporelle du mouvement brownien est markovienne, et donc entièrement déterminée par son propagateur, c’est-à-dire la (densité de) probabilité d’atteindre $x_1$ à la date $ t_1$ en partant de $ x_0$ à la date $ t_0$. Il est bien connu que ce propagateur est gaussien. Cependant, pour des applications pratiques, il est souvent nécessaire de restreindre le processus de Wiener à une demi-droite ou un intervalle, avec des conditions aux bords imposées (absorption ou réflexion). Dans l’exemple du prix d’un stock, ces dernières décriront les options d’achat et de vente. Dans cet article, je vais obtenir le propagateur du mouvement brownien avec dérive et conditions aux bords linéaires, avec une approche que j’espère pédagogique.

Propagateurs sans dérive

Considérons d’abord le mouvement brownien standard $ W(t)$ sans dérive. Sans perte de généralité, supposons $W(0)=0$. $W$ vérifie l’équation de Langevin :

\[\dot{W}(t) = \xi(t)\]

où $ \xi(t)$ est un bruit blanc gaussien avec des corrélations

\[\overline{\xi(t)\xi(t’)} = 2\sigma \delta(t-t’).\]

Avec ces conventions, le propagateur libre (c’est-à-dire sans conditions aux bords) $ P(x,t)$ est donné par la solution de l’équation de Fokker-Planck

\[ \partial_t P(x,t) = \sigma \partial_x^2 P(x,t)\]

avec la condition initiale $ P(x,0) = \delta(x)$. Cette équation aux dérivées partielles, également connue sous le nom d’équation de la chaleur, se résout facilement en prenant une transformée de Fourier. On trouve

\[ P(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}e^{-\frac{x^2}{4\sigma t}}.\]

Maintenant, intéressons-nous au propagateur avec une frontière absorbante en $ x=b > 0$. Dans l’équation de Fokker-Planck, ceci se traduit par $ P(b,t)=0$, ce qui rend plus difficile la transformation de Fourier. Cependant, on peut utiliser la méthode des images pour trouver la position : $P(b,t)=0$ est automatique en ajoutant une source négative en $ x=2b$ (symétrique de la source initiale par rapport à la frontière en $ b$). Cela revient à prendre la condition initiale $ P(x,0) = \delta(x)-\delta(x-2b)$. Finalement, le propagateur avec frontière absorbante en  $x=b$ est

\[P^{(b)}(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}\left[e^{-\frac{x^2}{4\sigma t}}-e^{-\frac{(x-2b)^2}{4\sigma t}}\right].\]

On traite le cas de deux frontières absorbantes en $ x=b>0$ et $ x=a<0$ de la même façon, mais il faut maintenant une infinité d’images, et le propagateur s’écrit sous forme d’une série qui peut être exprimée en termes des fonctions $\vartheta$ de Jacobi.

Propagateurs avec dérive

Généralisons maintenant au mouvement brownien avec une dérive $ \mu$. L’équation de Langevin pour $ W(t)$ devient

\[ \dot{W}(t) =\ mu + \xi(t).\]

Comme ci-dessus, le propagateur libre se déduit de l’équation de Fokker-Planck :

\[P(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}e^{-\frac{(x-\mu t)^2}{4\sigma t}}.\]

Introduisons maintenant une frontière absorbante en $ x=b$. La méthode des images n’est plus évidente à appliquer ! A cause de la dérive, le poids d’une trajectoire partant de $ x=0$ vers la frontière en $ x=b$ ne sera pas le même que celui de son image, de $ x=2b$ à $ x=b$. Cependant, dans le cas d’une dérive constante, une relation simple lie les poids avec et sans dérive. Il me semble que la façon la plus simple de s’en rendre compte est l’utilisation des intégrales de chemin. Le propagateur est donné par

\[ P(x_f,t) = \int_{x(0)=0}^{x(t)=x_f} \mathcal{D}[x]e^{-\int_0^t \mathrm{d}s, \frac{1}{4\sigma}\left(\dot{x}(s)-\mu\right)^2}\]

En développant maintenant “l’action” dans l’exponentielle, et en utilisant le fait que $ \mu$est constant, ainsi que les conditions aux bords, on trouve

\[ P(x_f,t) = e^{-\frac{\mu}{2\sigma}x_f+ \frac{\mu^2}{4\sigma}t}\int_{x(0)=0}^{x(t)=x} \mathcal{D}[x]e^{-\int_0^t \mathrm{d}s, \frac{1}{2\sigma}\left(\dot{x}(s)\right)^2}.\]

On obtient alors un poids qui ne dépend que de la position finale (hors de l’intégrale), et l’intégrale restante ne présente plus de dérive, ce qui signifie que nous connaissons déjà la solution, avec ou sans frontières ! Dans la littérature mathématique, vous trouverez fréquemment cette manipulation sous le nom de “Théorème de Cameron-Martin-Girsanov“, mais je trouve que l’approche via l’intégrale de chemin est beaucoup plus claire pour un physicien. Notons que si la dérive $ \mu$ dépend du temps, on ne peut pas sortir le poids de l’intégrale de chemin, car il dépend de toute la trajectoire, et pas seulement du point final. C’est pourquoi il s’agit d’un problème bien plus difficile, bien que le propagateur libre soit toujours trivial à écrire !

La formule finale pour le propagateur du mouvement brownien avec dérive $ \mu$ et une frontière absorbante en $ x=b$ (connu sous le nom de “formule de Bachelier-Levy”) est donc

\[P^{(b,mu)}(x,t) = \frac{1}{\sqrt{4\pi\sigma t}}e^{-\frac{mu}{2\sigma}x+ \frac{mu^2}{4\sigma}t}\left[e^{-\frac{x^2}{4\sigma t}}-e^{-\frac{(x-2b)^2}{4\sigma t}}\right]. \]

C’est tout ce dont nous avons besoin pour calculer des choses comme des temps de premier passage, des probabilités de survie, etc. La généralisation à deux frontières absorbantes peut être déduite de la solution sans dérive, en multipliant par le même poids.

J’espère vous avoir fait sentir qu’avec les bons outils, obtenir ces propagateurs n’a rien de miraculeux. Si vous avez des questions ou des commentaires, je serai ravi de les entendre ! Pour approfondir, diverses voies sont envisageables, comme des généralisations des méthodes exposées ici aux processus de Ornstein-Uhlenbeck, de Bessel, ou encore des conditions aux bords plus compliquées.

Merci de m’avoir lu !

Epidemias, parte I

Hardcore

Trombei esses dias com um artigo de 2001 sobre propagação de vírus em humanos e sua comparação a vírus de internet. O estudo da disseminação de infecções é velho na física, diversos modelos de epidemias são estudados e muitos com resultados razoavelmente bons. Um ponto comum em modelos estudados pela epidemiologia é a noção de threshold, ou seja, o valor mínimo de eficiência que uma infecção deve ter para conseguir se propagar. Novamente, esse é um post hardcore, os corajosos sigam-me, vamos ver com cuidado como isso funciona.

Esse post ficou grande, então parti em dois. Nesse primeiro, comento o modelo de propagação de vírus em seres humanos, no próximo trato da internet.

Representamos indivíduos como pontos e o contato entre indivíduos como linhas em um grande grafo. Passamos a nos perguntar qual seria um modelo eficaz de grafo para modelizar as relações humanas. Certamente não conhecemos pessoas aleatoriamente na rua, costumamos conhecer melhor nossos vizinhos, familiares e colegas de trabalho, que por sua vez conhecem-se. Um jeito legal de modelizar essa ideia, que chamamos de grafo de small world, é tomar um anel de pontos, sendo cada ponto ligado a seus $ m$ vizinhos mais próximos. Em seguida, passamos “em revista” essas conexões da seguinte forma: cada ponto terá uma probabilidade $ p$ de desfazer uma conexão com um vizinho seu no sentido horário e associá-la a algum outro ponto aleatoriamente. A ideia do sentido horário é apenas para evitar mexer na mesma conexão duas vezes. Se $ p=1$, essa revista resulta em um grafo em forma de anel com conexões completamente aleatórias. Se $ p=0$, o grafo continuará no modelo do conhecimento apenas da vizinhança, mas, se $ p$ é pequeno, o grafo tornar-se-á um modelo muito próximo do que queremos: vizinhos possuem uma chance maior de se conhecerem do que conhecerem alguém do outro lado do grafo, mas de vez em quando algum vizinho seu conhece alguém lá longe, o que é bem coerente com a realidade e até permite aquelas coincidências do tipo “hoje meu avô almoçou com o pai de um professor meu”.

O resultado está bem representado nessa figura, onde os pontos são inicialmente ligados a seus $ 2m$ vizinhos mais próximos, com $ m=2$. É o famoso modelo de Watts-Strogatz, ou modelo de small-world:

Podemos, com esse modelo bonitinho, estudar a propagação de infecções entre seres humanos. Imaginemos um vírus que, a cada contato entre são e infectado, tenha uma probabilidade $ nu$ de contágio e, a cada “turno”, tenha uma probabilidade $ \delta$ de cura de um infectado. A chance de um ponto possuir $ k$ vínculos é $ P(k)= \frac{1}{\langle k\rangle}e^{-k/\langle k\rangle}$ nesse modelo, e $ \langle k\rangle=2m$. A equação estocástica dessa criança não é bonita, então vamos usar um mean field, calcular o que vai acontecer com a densidade média de infectados supondo que todos os pontos possuem o mesmo número de conexões, $ 2m$, naturalmente, e que todos estão ligados a um número de infectados igual ao valor médio de infectados, essas são as características de uma abordagem de campo médio. A densidade total de infectados, denotada $ \rho_t$ por ser uma função do tempo, é apenas o número total de infectados dividido pelo número total de pessoas. Diremos que a variação na densidade de infectados $ \rho_t$ (a fração da população infectada naquela “rodada”), que toma valores entre 0 e 1, perde a cada “rodada” o equivalente a $ \delta\rho_t$ (os que se curaram) e ganha o equivalente ao contato são-infectado, que é $ nu\langle k\rangle\rho_t(1-\rho_t)$, ou seja, chance de um contato ser são-infectado $ \rho_t(1-\rho_t)$ multiplicada pelo número médio de contatos $ \langle k\rangle$ e pela chance de esse contato transmitir a infecção $ \nu$. Dividindo a equação toda por $ \delta$ e sem precisar ir muito além de uma reparametrização no tempo, podemos escrever:

\[\frac{d}{dt}\rho_t=-\rho_t+\lambda\langle k\rangle\rho_t(1-\rho_t).\]

Com $ \lambda= \frac{\nu}{\delta}$ o valor que interessa para medir a eficácia de uma infecção. E essa equação, ainda que seja apenas uma teoria de campo médio, já nos dá bastante informação sobre o que está acontecendo. Você até pode pedir para o Wolfram Alpha resolver essa para você, eu não faria diferente, mas vale mais estudar algumas propriedades na mão. Chamamos densidade de persistência o valor de $ rho_t$ para o qual sua derivada temporal é zero. Se ela não é nula, temos uma infecção persistente em uma fração da população. Igualando a equação acima a zero, percebemos duas respostas possíveis: ou $ rho_t = 0$, ou $ \rho_t=1- \frac{1}{\lambda\langle k\rangle}$. Como a densidade de infectados não pode ser menor que 0, naturalmente, descobrimos que, se $ \lambda < \frac{1}{\langle k\rangle}$, então a única solução possível é $ \rho_t=0$, pois a segunda solução fica absurda. Mas, se a infecção é mais eficaz que o threshold $ \frac{1}{\langle k\rangle}$, teremos um crescimento da densidade de persistência que tende a 1 quando $ \lambda$ se aproxima do infinito.

Densidade de persistência com m=2.

E fica claro que o threshold, dependendo do inverso do número médio de contatos, que nesse modelo é proporcional ao número de vizinhos, nos dá uma informação esperada: quanto maior o número de vizinhos (contatos), menos eficaz uma doença deve ser para se propagar. Se uma doença é bem eficaz, podemos reduzir seu impacto com uma diminuição nos vizinhos, se está com gripe, não saia de casa.

Esse modelo é bem coerente com diversas observações de epidemias simples em populações. Por esse mecanismo de small-world, apenas doenças suficientemente poderosas, com eficácia $ \lambda$ maior que um determinado limite, podem se propagar pela população, possuem vida no longo prazo. A mudança de fase ocorre no threshold $ \lambda^\star= \frac{1}{\langle k\rangle}$ que, no gráfico acima, ocorre em $ \lambda = 0,25$. O próximo passo nesse modelo seria estudar um grafo em que pessoas podem morrer, nascer, tornar-se imunes, tomarem vacinas ou tornarem-se mais suscetíveis a doenças com o tempo.

Tudo isso é muito interessante, mas não consegue explicar a propagação de vírus na internet. Em um próximo post, veremos que, na internet, small-world não pode ser aplicado (afinal, todos conhecem o Google ou o Facebook) e a própria noção de threshold, por um fenômeno estatístico fascinante, desaparece.

Prova surpresa

Geek Hardcore Rookie

O post de hoje não é muito informativo, é mais perturbador, um probleminha de lógica para não deixar você dormir.

Um professor chega à aula e avisa sua classe que aplicará, no mês de abril, uma prova surpresa. Um aluno, que futuramente se tornaria matemático na área de lógica, pergunta o que o professor entende por surpresa. Estranhando a pergunta, o professor responde que uma prova surpresa é uma prova tal que, no início de cada dia do mês de abril, os alunos não podem deduzir ou saber que a prova será aplicada naquele dia. Triunfante, o aluno conclui que não haverá prova nenhuma, pois:

– Se a prova não for aplicada até o 30 de abril, então ela terá que necessariamente ser nesse dia e, no começo do dia, já não será surpresa. Então o 30 de abril não pode ser a data da prova.

– Se a prova não for aplicada até o 29 de abril, então ela terá que ser necessariamente nesse dia, pois o 30 está excluído. Mas, assim sendo, ela não será uma surpresa, então o dia 29 não pode ser o dia da prova.

– Se a prova não for aplicada até o 28 de abril, então ela terá que ser necessariamente nesse dia, pois o 30 e o 29 estão excluídos. Mas, assim sendo, ela não será uma surpresa, então o dia 28 não pode ser o dia da prova. E assim por diante, ele exclui todos os dias.

No dia 12, ele recebe a prova. Onde ele errou?

O gás de Coulomb-Dyson

Hardcore

Esses dias trombei com um assunto bonito, uma parte do estudo de matrizes aleatórias que realmente achei interessante. Claro, trabalho com matrizes aleatórias, então é um pouco normal eu encontrar coisas bonitas aqui e ali, mas dificilmente algo tão bonito quanto o que achei outro dia em um desses livros de teoria espectral. Aviso, esse post é nível hardcore e provavelmente o mais intenso que já postei até agora, recomendo discrição. Se você não é um físico ou matemático, não vai pescar muita coisa do texto, não aconselho sua leitura.

As matrizes aleatórias servem para bastante coisa. Diversos modelos envolvendo um operador com perturbação aleatória e fora de nosso controle podem ser colocado em um formato de matriz aleatória, a transmissão de dados de um conjunto de antenas a outro com ruído pode ser modelizado por uma equação linear com uma matriz aleatória também, muita coisa entra nessa categoria e em cada um desses estudos a pergunta é recorrente: jogando entradas aleatórias em uma matriz, com uma densidade de probabilidade $ P(x)$ para as entradas, qual será a densidade de probabilidade de seus autovalores?

Porque, no fundo, é nisso que estamos interessados. Na quântica eles serão as medidas possíveis, nas equações diferenciais lineares eles nos darão a estabilidade do sistema, autovalores são a alma das matrizes. Mas se encontrar os autovalores de uma matriz $ n \times n$ bem definidas já não é tarefa tão fácil, resolver polinômios sempre dá preguiça, é difícil não sofrer só ao pensar como será com matrizes cuja única informação é a probabilidade de obter suas entradas entre dois valores.

Primeiro vamos estabelecer o que eu quero dizer com “probabilidade de uma matriz”. De forma civilizada, eu precisaria temperar esse texto com alguns detalhes da medida usada, do espaço em questão, mas isso é um blog e não um artigo; essa história está mais bem contada neste excelente artigo sobre matrizes aleatórias, ainda nas primeiras páginas. Tomemos um exemplo fácil, um vetor $ (x,y)$. Falar de sua densidade de probabilidade, é falar de uma função $ P(x,y)$ tal que, para descobrir a chance de encontrar esse vetor com as coordenadas $ 0\leq x \leq 3$ e $ 0 \leq y \leq 2$ será o resultado da integral $ \int_0^2 dy \int_0^3dx P(x,y)$, simples assim. Se $ x$ e $ y$ são independentes, certamente teremos $ P(x,y)=P(x)P(y)$, mas isso não será necessariamente verdade. Falar de densidade de probabilidade de uma matriz é falar da densidade conjunta de suas entradas, ou seja, sua j.p.d.f. (joint probability density function). Eu usaria o termo em português, mas não gosto de abreviar função densidade de probabilidade.

É claro que nem as matrizes nem as probabilidades podem ser quaisquer, problemas gerais demais não levam a lugar nenhum. Vou me restringir ao caso real, tudo pode ser generalizado para complexo com as devidas trocas. Listo dois tipos de matrizes cujo estudo das probabilidades é interessante: matrizes cujas probabilidades das entradas são independentes e matrizes que possuem a probabilidade invariante por conjugação.

Essa última propriedade é mais sutil, mas é simples. Dizer que a j.p.d.f. de uma matriz $ M$ é invariante por conjugação é dizer que, para toda $ U$ não singular, teremos que $ P(M)=P(UMU^{-1})=P(M’)$, ou seja, a chance de obter um elemento da classe de conjugação de $ M$ é a mesma chance de se obter $ M$.

E essa propriedade serve para uma manobra bem útil. Se a j.p.d.f. é invariante por conjugação, e se eu consigo diagonalizar a matriz, a probabilidade da matriz será a mesma probabilidade de sua forma diagonal, com seus autovalores como entradas, o que me permite de maneira fácil obter a distribuição de probabilidade de seus autovalores. Por vários motivos, físicos e matemáticos, gostamos de estudar matrizes simétricas, ou hermitianas se complexas. Isso vai garantir a diagonalização por matrizes unitárias e a existência de um alegre conjunto de autovalores bem reais.

Teorema: O único grupo de matrizes aleatórias invariantes por conjugação unitária e cujas entradas possuem p.d.f. independentes é o grupo das matrizes gaussianas, cujas entradas possuem como p.d.f. uma distribuição normal.

Nem arrisco tentar demonstrar isso, tomaria este post e mais outros três. Esse teorema nos inspira a aprofundar nosso estudo das matrizes gaussianas. Como elas são invariantes por conjugação unitária (estamos ainda com matrizes simétricas se reais e hermitianas se complexas, então elas são diagonalizáveis por uma matriz unitária), podemos escrever que $ P(M) = P(D)$, onde $ D$ é sua forma diagonalizada. Para que isso seja possível, vamos nos restringir às matrizes gaussianas simétricas (hermitianas se são complexas). Como em $ P(M)$ as probabilidades das entradas individuais são independentes, nosso instinto nos diz que em $ P(D)$ as coisas serão parecidas, que os autovalores terão probabilidades independentes, e não poderíamos estar mais errados. Porque escrever $ P(D)$ é escrever $ P(\lambda_1,\lambda_2,ldots,\lambda_n)$, uma probabilidade que depende dos autovalores, isso é uma mudança de variável e, como a probabilidade sempre se dá integrando essa densidade de probabilidade, precisamos levar em conta o jacobiano dessa transformação que, para nosso desespero, acopla todos os autovalores. Tal probabilidade já é conhecida há algum tempo, a j.p.d.f. dos autovalores da matriz gaussiana:

\[P(\lambda_1,ldots,\lambda_n) = C_k e^{-\beta \sum_k \lambda_k^2} \prod_{j<k}|\lambda_k-\lambda_j|^\beta \]

E o último termo da direita, a parte do jacobiano relativa aos autovalores, não é ninguém menos que o determinante de Vandemonde. O $ \beta$ é um valor referente ao tipo da matriz gaussiana, vale 1 se ela é real, 2 se é complexa e 4 se estamos nos quatérnions. E agora vem a parte bonita: Dyson (nisso fui corrigido, disseram ser Wigner, deixo a polêmica) percebeu que essa j.p.d.f. poderia ser colocada de uma forma mais familiar, bastava apenas jogarmos o determinante de Vandermonde para o expoente como $ \prod_{j<k}|\lambda_k-\lambda_j|^\beta = e^{\beta \sum_{j<k} \log |\lambda_k-\lambda_j|}$, teremos que a probabilidade será um múltiplo de uma grande exponencial. Chamar aquele número de $ \beta$ é extremamente sugestivo. Um leitor atento já deve ter percebido, contemplamos um peso de Boltzmann, em uma analogia perfeita a uma j.p.d.f. do sistema canônico:

\[P(\lambda_1,ldots,\lambda_n) = C_k e^{-\beta\left(\sum_k\lambda_k^2-sum_{j<k}\log |\lambda_k-\lambda_j|\right)} = \frac{1}{Z}e^{-\beta H}.\]

A magia dessa interpretação é poder importar todas as ferramentas da física estatística para resolver esse intrincado problema de álgebra linear. Se imaginarmos que cada autovalor representa a posição de um elétron confinado a uma linha (que representará o eixo real), atraídos ao centro por um potencial harmônico (o termo em $ \sum_k\lambda_k^2$ ) e submetidos à repulsão coulombiana mútua (que em sua forma bidimensional é $ \log |x_i-x_j|$), teremos um sistema físico cujas posições dos elétrons são equivalentes às posições dos autovalores da matriz gaussiana. E, caramba!, isso é muito bonito.

A analogia não é apenas formal, podemos extrair diversas propriedades dessa distribuição com técnicas do ensemble canônico (meu orientador, aliás, fez a carreira dele nisso). E este é um de meus exemplos favoritos de física ajudando matemática, ainda que, se fôssemos contar pontos nisso, a competição seria injusta, estaríamos perdendo por uma boa margem.

Introdução

Geek Hardcore Rookie

Primeiro post, vamos ver no que isso vai dar.

Todo texto ou artigo começa com uma introdução, aquela parte que quase ninguém lê. Precisava ser justo e começar eu também nessa, explicando a que veio este blog.

Gosto de física e de matemática. Sou físico estatístico, terminando meu mestrado e caminhando para um doutorado, se os ventos permitirem. Muita coisa me fascina no que faço, sou desses apaixonados pela ciência, pelos resultados, pelo desafio e pelas pessoas que a fazem. Gosto de tanta coisa que esse blog ficaria injusto se focasse em um público, se fosse apenas divulgação ou apenas física avançada, então decidi criar três categorias e postar de acordo. A primeira, rookie, é para quem não viu mais ciência do que lhe foi martelado no colegial com musiquinhas e provas de raciocínio a lápis e resposta a caneta. A segunda, geek, é para os que se aventuraram a continuar, ou ainda estão neste treinamento, nas artes das exatas pelo ensino superior; penso em engenheiros, químicos ou estudantes de faculdade nessa categoria. A última, hardcore, é para os que concluíram o treinamento em física ou matemática, estão à vontade com epsilons e deltas, com bras e kets, com difeomorfismos e fibrados tangentes.

Não tenho tema específico, gosto de coisas bonitas na ciência. Sou físico estatístico, o que explica o título do blog; nessa área, estamos cansados de fazer o que chamamos de “soma sobre todas as configurações possíveis”, a base da definição de entropia e de algumas coisinhas mais. Claro, os posts podem acabar tendenciosos, mas não posso evitar, não posso escrever sobre super-cordas se não aprendi super-cordas, não me arrisco nas álgebras de Hopf se jamais me aventurei tanto na teoria das representações. Talvez mude de ideia, talvez Hopf me conquiste, não tenho compromissos com áreas ou gostos no que escrevo.

Vejamos no que estes textos vão dar. Meu único propósito é divertir quem lê, vocês, tentando compartilhar um trecho da beleza na ciência, que tanto é profanada nos impiedosos bloquinhos em plano inclinado de nosso ensino médio. Porque toda ciência começa como filosofia e termina como arte, sua beleza é difícil de deixar de perceber. Ela é, diferentemente de toda poesia, verdade, o que a reveste de um ar selvagem, belo e a ser desvendado. Mais bela que poemas, a ciência é a poesia da realidade, escrita no que somos e vemos, traduzida em matemática, um jogo de xadrez com o divino, em que temos como única missão entender as regras e jogar.