# Minuit, l’heure du crime ? heure d’hiver ou d’été ? Suite à mon précédant billet où j’essayais d’étudier la dynamique des appels aux 911 à Montréal, j’avais noté qu’un seul type d’appels avait lieu aux alentours de minuit, à savoir les trouble à la paix. Mais minuit est une notion quelque peu relative. En tous les cas, c’est la réflexion que je me faisais alors que les jours commencent à rallonger… car l’hiver à Montréal, minuit ça peut presque 8 heures après le coucher du soleil, et l’été à peine 4 heures…C’est aussi ce que me demandait @simon dans un commentaire.

En fait, de manière assez étrange, les cycles sont relativement stables en fonction des saisons (avec dans le cercle au centre l’ensoleillement de la journée – avec en bleu la nuit) Je mettrais peut être cela sur le compte que Montréal est relativement proche de l’heure solaire, au sens où midi est assez souvent le milieu de la journée et que minuit est le milieu de la nuit.

Dans le même état d’esprit, on peut se demander s’il y a un impact du changement d’heure. Par exemple si on considère le changement d’heure qui survient au printemps. Le graphique de gauche représente les appels le week-end qui précédait le changement d’heure (vendredi, samedi et dimanche qui précédaient) et à droite le week-end qui a suivi (même trois jours),  Il y a peut être un léger désaxage vers la droite, mais à peine…

Et si on regarde d’autres types d’appels, comme les coups de feu, on observe le début de la période d’appel (un peu plus intense) suit le commencement de la nuit: s’ils finissent vers 4 heures du matin, les appels arrivent plus tôt l’hiver que l’été en revanche, on a un impact presque nul du changement d’heure,  On peut aussi regarder les conflits, qui présentaient une dynamique relativement différente, on observe ici une relative stabilité et encore une fois, on a un impact presque nul du changement d’heure,  # ma que calor !

J’avais commencé plusieurs analyse (ici ou ) sur les appels passé au 911, à la police, à Montréal, en fonction de l’heure et des jours de la semaine afin de détecter d’éventuels cycles. Mais ces derniers jours, maintenant qu’il commence à faire beau et chaud, je me suis demandé s’il n’était pas possible de croiser avec la température… Car j’aurais l’idée naïve qu’on ne va pas commettre un cambriolage quand il fait -25°C, ou que les conflits sont en hausse s’il fait +35°C… Si on regarde assez rapidement, en comptant le nombre d’appels en fonction de la température, et en faisant une régression splines, on obtient autrement dit, entre -20°C et +30°C, il y a environ 50% d’appels en plus. Ce qui est loin d’être négligeable ! Mais on peut aussi croiser la température avec l’heure d’appel, avec l’heure en abscisse et la température en ordonnée  Aussi, la région ou le nombre d’appels est le plus élevé est indiquée ci-dessous: on note que l’on a plus d’appels quand il faut chaud, et plutôt en fin de journée, Pour les conflits, par exemple, on a davantage de conflits en fin d’après midi comme on l’avait relevé, mais surtout, ces derniers sont très liés à la température: il y a d’autant plus d’appels qu’il faut chaud.  Les appels sont d’autant plus importants dans la région bleue ci-dessous, Si l’on se concentre uniquement sur la température et les conflits, on a l’évolution suivante après presque deux fois plus d’appels s’il fait +30°C par rapport à -20°C. De même pour les agressions sexuelles, qui augmentent très fortement avec la température (environ +25%), En revanche les holds ups sont indépendants de la température, En revanche les alarmes cambrioleurs semblent générer plus d’appels s’il fait chaud, mais aussi s’il fait froid… # 911, jour après jour

Après deux billets (ici puis ) sur les cycles intrajournaliers des appels au 911, on peut se demander comment les crimes évoluent au cours de la semaine. Pour l’ensemble des appels passés au 911, on a la distribution suivante i.e. un pic les vendredi soir et samedi soir, et un creux le dimanche. Si on regarde les appels pour des cambriolages, on a la distribution suivante avec des pics en matinée, les vendredi après midi, et les fins de semaine. On peut aussi suivre les troubles de la paix, qui surviennent certes vers minuit, mais essentiellement en fin de semaine. Ce qui contraste assez avec les hold-ups, Manifestement, il y a des tendances assez claires. La prochaine étape sera de regarder un peu les saisons, ou mieux, l’impact du climat… à suivre donc…

# Circular or spherical data, and density estimation

I few years ago, while I was working on kernel based density estimation on compact support distribution (like copulas) I went through a series of papers on circular distributions. By that time, I thought it was something for mathematicians working on weird spaces…. but during the past weeks, I saw several potential applications of those estimators.

• circular data density estimation

Consider the density of an angle say, i.e. a function such that with a circular relationship, i.e. . It can be seen as an invariance by rotation.
von Mises proposed a parametric model in 1918 (see here or there), assuming that where is Bessel modified function of order 1, (which is simply a normalization parameter). There are two parameters here, (some concentration parameter) and mu a direction.
From a series of observed angles , the maximum likelihood estimator for kappa is solution of where and and where , where those functions are modified Bessel functions. Well, that estimator is biased, but it is possible to improve it (see here or there). This can be done easily in R (actually Jeff Gill – here – used that package in several applications). But I am not a big fan of that technique….

• density estimation for hours on simulated data

A nice application can be on the estimation of the daily density of a temporal events (e.g. phone calls as we’ll see later on, or email arrival time). Let is the time (in hours) for the th observation (the th phone call received). Then set The time is now seen as an angle. It is possible to consider the equivalent of an histogram,

```set.seed(1)
library(circular)
X=rbeta(100,shape1=2,shape2=4)*24
Omega=2*pi*X/24
Omegat=2*pi*trunc(X)/24
plot(Ht, stack=FALSE, shrink=1.3, cex=1.03,
points(Ht, rotation = "clock", zero =c(rad(90)),
col = "1", cex=1.03, stack=TRUE )

rose.diag(Ht-pi/2,bins=24,shrink=0.33,xlim=c(-2,2),ylim=c(-2,2),
axes=FALSE,prop=1.5)```  or a kernel based estimation of the density (the gray line on the right).

```circ.dens = density(Ht+3*pi/2,bw=20)
plot(Ht, stack=TRUE, shrink=.35, cex=0, sep=0.0,
axes=FALSE,tol=.8,zero=c(0),bins=24,
xlim=c(-2,2),ylim=c(-2,2), ticks=TRUE, tcl=.075)
lines(circ.dens, col="darkgrey", lwd=3)
text(0,0.8,"24", cex=2); text(0,-0.8,"12",cex=2);
text(0.8,0,"6",cex=2); text(-0.8,0,"18",cex=2)```

The c ode looks rather simple. But I am not very comfortable using codes that I do not completely understand. So I did my own. The first step was to get a graph similar to the one we have on the right, except that I prefer my own kernel based estimator. The idea is that instead of estimating the density on , we estimate it on the sample . Then we multiply by 3 to get the density only on . For the bandwidth, I took the same as the one that we would have taken on The code is simply the following

```U=seq(0,1,by=1/250)
O=U*2*pi
U12=seq(0,1,by=1/24)
O12=U12*2*pi
X=rbeta(100,shape1=2,shape2=4)*24
OM=2*pi*X/24
XL=c(X-24,X,X+24)
d=density(X)
d=density(XL,bw=d\$bw,n=1500)
I=which((d\$x>=6)&(d\$x<=30))
Od=d\$x[I]/24*2*pi-pi/2
Dd=d\$y[I]/max(d\$y)+1

plot(cos(O),-sin(O),xlim=c(-2,2),ylim=c(-2,2), type="l",axes=FALSE,xlab="",ylab="") for(i in pi/12*(0:12)){ abline(a=0,b=tan(i),lty=1,col="light yellow")} segments(.9*cos(O12),.9*sin(O12),1.1*cos(O12),1.1*sin(O12)) lines(Dd*cos(Od),-Dd*sin(Od),col="red",lwd=1.5) text(.7,0,"6"); text(-.7,0,"18") text(0,-.7,"12"); text(0,.7,"24") R=1/24/max(d\$y)/3+1 lines(R*cos(O),R*sin(O),lty=2)```

Note that it is possible to stress more (visually) on hours having few phone calls, or a lot (compared with an homogeneous Poisson process), e.g.

```plot(cos(O),-sin(O),xlim=c(-2,2),ylim=c(-2,2),
type="l",axes=FALSE,xlab="",ylab="")
for(i in pi/12*(0:12)){
abline(a=0,b=tan(i),lty=1,col="light yellow")}
segments(2*cos(O12),2*sin(O12),1.1*cos(O12),1.1*sin(O12), col="light grey")
segments(.9*cos(O12),.9*sin(O12),1.1*cos(O12),1.1*sin(O12))
text(.7,0,"6")
text(-.7,0,"18")
text(0,-.7,"12")
text(0,.7,"24")
R=1/24/max(d\$y)/3+1
lines(R*cos(O),R*sin(O),lty=2)
AX=R*cos(Od);AY=-R*sin(Od)
BX=Dd*cos(Od);BY=-Dd*sin(Od)
COUL=rep("blue",length(AX))
COUL[R<Dd]="red"
CM=cm.colors(200)
a=trunc(100*Dd/R)
COUL=CM[a]
segments(AX,AY,BX,BY,col=COUL,lwd=2)
lines(Dd*cos(Od),-Dd*sin(Od),lwd=2)```

We get here those two graphs,  To be honest, I do not really like that representation – even if it looks nice. If we compare that circular representation to a more classical one (from 0:00 till 23:59 one the graph on the left, below), I do have a problem to interpret the areas in blue and pink.  density of wind direction

On the left, we compare two densities, so the area in pink is the same as the area in blue. But here, it is no longer the case: the area in pink is always larger to the one in blue. So it might help so see when we have a difference, but there is a scaling issue that we cannot discuss further… But less us see if we can use that estimation technique to several problems. A standard application when studying angles is wind direction. For instance, in Montréal, it is possible to find hourly observations, starting in 1974 (we just need a R robot to pick up the information, but I’ll tell more about that in another post, someday). Here, we have directly an angle. So we can use a code rather similar to the one used above to estimate the distribution of wind direction in Montréal.  density of 911 phone calls

Note that our estimate is consistent with several graphs that can be found on meteorological websites (e.g. the one above on the right, that was found here).

In a recent post (here) I wanted to check about the “midnight crime” myth, using hours of 911 phone calls in Montréal. That was for all phone calls. But if we look more specifically, for burglaries, we have the distribution on the left, and for conflicts the one on the right  We do clearly observe that gun shots occur a bit before midnight. See also here for another study, but this time in NYC (thanks @PAC for the link).while for gun shots, we have the distribution on the left, and for “troubles” (basically people making too much noisy in parties) or “noise” the one on the right

• density of earth temperatures, or earthquakes

Of course it is also possible to work in higher dimension. Before, we went from densities on to densities on the unit circle . But similarly, it is possible to go from to the unit sphere . A nice application being global climate studies, The idea being that point on the left above are extremely close to the one on the right. An application can be e.g. on earthquakes occurrence. Data can be found here.

```library(ks)
X=cbind(EQ\$Longitude,EQ\$Latitude)
Hpi1 = Hpi(x = X)
DX=kde(x = X, H = Hpi1)
library(maps)
map("world")
points(X,cex=.2,col="blue")
Y=rbind(cbind(X[,1],X[,2]),cbind(X[,1]+360,X[,2]),
cbind(X[,1]-360,X[,2]),cbind(X[,1],X[,2]+180),
cbind(X[,1]+360,X[,2]+180),cbind(X[,1]-360,X[,2]+180), cbind(X[,1],X[,2]-180),cbind(X[,1]+360, X[,2]-180),cbind(X[,1]-360,X[,2]-180)) DY=kde(x = Y, H = Hpi1) library(maps) plot (DY,add=TRUE,col="purple")```

Without any correction, we get the red level curves. The pink one integrates correction. # Minuit, l’heure du crime ?

Il y a quelques jours, j’étais invité par Rémi pour me rendre au poste de police pour discuter de la possibilité de consulter quelques données. Et Rémi m’a donné un historique de près d’un million d’appels téléphonique d’urgence, sur Montréal, passés au cours des six dernières années. Histoire de prendre un peu en main les données, j’ai voulu regarder (inspiré par une remarque de Rémi) si minuit était vraiment une heure de crime1

Ce dont je dispose ce sont les heures d’appels d’urgence (le fameux 911), avec un code (formellement, le code indiqué par la personne lors de l’appel, pas celui établi par la patrouille une fois sur place). Globalement, la distribution des appels (par seconde) au cours de la journée ressemble à ça, avec les six années en noir, 2005 en bleu et 2010 en rouge (histoire de visualiser une stabilité temporelle).
On peut aussi se focaliser sur les causes données lors de l’appel, afin de voir quel criminel commet un crime à minuit. Par exemple, pour les alarmes consécutives à un cambriolage, on obtient la répartition suivante, Manifestement, minuit est un creux en matière de cambriolage. N’en déplaise à Arsène Lupin. Les hold-ups, eux aussi surviennent pendant la journée, Les intrusions-effractions aussi Sinon, de manière générale, les infractions au code criminel surviennent en fin de journée (mais baissent en fin de soirée) et les conflits, suivent la même distribution Autrement dit, les conflits surviennent quand on rentre du travail… mais pas à minuit. En revanche, si on commence à regarder les appels pour signaler descoups de feux, ils surviennent un peu après les conflits… mais un peu avant minuit. Ah, mais finalement je crois que j’ai trouvé le crime qui survient à minuit, Et dans ma base, cela correspond au code 063, ce qui signifie trouble à la paix. Autrement dit, si un appel arrive au 911 à minuit, c’est rarement pour signaler un (vrai) crime (au sens de ceux qu’on trouve dans les polars), mais plus pour signaler un problème de beuverie (que l’on peut aussi trouver chez Ian Rankin d’ailleurs). Damned ! Moi qui pensait résoudre des crimes avec mes séries de chiffres…
.
1 Je n’ai pas réussi à trouver d’où pouvait sortir cette idée reçue. En lisant le dernier Ian Rankin (ou plutôt la dernière enquête de John Rebus), je me suis fait la réflexion que dans beaucoup d’enquêtes, les policiers arrivaient vers minuit. J’ai eu l’impression que c’était aussi souvent le cas dans bon nombre d’épisodes de séries télé “policières” (CSI, Castle, Bones, ). En fait si quelqu’un sait où je pourrais trouver une base avec des dates de décès des personnages dans les séries policières ou des polars, je suis preneur.