FR 2967512 A1 20120518 FR 1059452 A 20101117 Procédé d'analyse de données longitudinales, programme d'ordinateur et système correspondants La présente invention concerne un procédé d'analyse de données longitudinales caractérisant l'évolution d'au moins une première variable en fonction d'au moins une deuxième variable, comprenant des étapes de détermination de sous-intervalles adjacents de variation d'au moins une desdites première et/ou deuxième variables et de caractérisation desdites données sur lesdits sous-intervalles. Elle s'applique en particulier à la réalisation automatique de graphiques d'évaluation de modèles mathématiques, tels que les graphiques VPC (pour « Visual Predictive Check »). Ce type de graphiques permet de comparer des données obtenues par simulation, au moyen d'un modèle mathématique, à des données réelles obtenues par des observations, lorsque ces données se présentent sous la forme de données longitudinales, c'est-à-dire expriment l'évolution d'une première variable, notée par la suite y;, en fonction d'une deuxième variable, notée par la suite x;. Une telle comparaison permet alors de choisir, parmi plusieurs modèles candidats, celui présentant les meilleures capacités de restitution des données et/ou de prédiction. Les graphiques VPC sont par exemple utilisés pour l'évaluation de modèles pharmacocinétiques/pharmacodynamiques (PK/PD), modélisant les évolutions respectives de la concentration d'une substance active dans un organisme et d'un effet de cette substance active sur l'organisme en fonction du temps. Une première étape de réalisation d'un graphique VPC consiste à générer un grand nombre de données simulées au moyen du modèle considéré, c'est-à-dire à calculer les valeurs prises par la première variable, ou grandeur étudiée, pour différentes valeurs de la deuxième variable, appelée par la suite coordonnée longitudinale, puis à représenter sur un même graphique la distribution statistique des données observées et celle des données simulées, dans le but de comparer ces distributions. En effet, une différence importante entre ces deux distributions traduit généralement une mauvaise adéquation du modèle au phénomène étudié. Les données observées et simulées sont ainsi représentées sur un graphique portant en abscisse la coordonnée longitudinale et en ordonnée les valeurs prises par la grandeur étudiée, sous la forme de deux ensembles de points. Ces deux ensembles ne peuvent généralement pas être comparés directement l'un à l'autre. En effet, les coordonnées longitudinales des données simulées peuvent différer de celles des données observées, et la grandeur étudiée peut prendre de multiples valeurs pour une même coordonnée longitudinale. 2 La comparaison des données observées et simulées est donc généralement réalisée en découpant l'intervalle de variation de la coordonnée longitudinale en une pluralité de sous-intervalles adjacents de variation, et en caractérisant les valeurs réelles et/ou simulées prises par la grandeur étudiée sur chacun de ces sous-intervalles, par exemple en calculant les percentiles de la grandeur étudiée sur chacun des sous-intervalles. La comparaison des percentiles des données observées et simulées permet alors de mesurer l'adéquation du modèle aux données réelles. Les résultats de cette comparaison dépendent du choix des sous-intervalles. Il est donc crucial de choisir des sous-intervalles permettant une caractérisation pertinente de la distribution des données. En particulier, ces sous-intervalles doivent être suffisamment larges, c'est-à-dire contenir suffisamment de données, pour que les valeurs réelles et/ou simulées prises par la grandeur étudiée sur ces sous-intervalles soient statistiquement significatives. Cependant, pour être en mesure de caractériser correctement l'évolution de la grandeur étudiée en fonction de la coordonnée longitudinale, il est nécessaire de déterminer un nombre suffisant de sous-intervalles, donc d'en limiter la largeur. Pour choisir ces sous-intervalles, une première méthode connue consiste à fixer arbitrairement le nombre Kde sous-intervalles, et à découper l'intervalle de variation de la coordonnée longitudinale en Ksous-intervalles de mêmes largeurs. Cette approche n'est pas satisfaisante, car les sous-intervalles obtenus ne dépendent pas de la distribution des données, distribution qui n'est généralement pas homogène, de telle sorte que certains sous-intervalles peuvent contenir un grand nombre de données tandis que d'autres sous-intervalles peuvent être vides. Pour équilibrer la répartition des données dans les différents sous-intervalles, une deuxième méthode consiste à déterminer des sous-intervalles présentant tous le même effectif. Cette méthode est plus satisfaisante que la précédente, car elle offre une répartition plus équilibrée des données dans les sous-intervalles. Cependant, cette méthode ne permet pas de contrôler la dispersion des données dans chacun des sous-intervalles. Un même sous-intervalle peut par exemple contenir des données de coordonnées x; très différentes, tandis que des données de coordonnées x; très proches peuvent être réparties dans des sous-intervalles différents. De plus, alors que la répartition d'un nombre N de données en K sous-intervalles de mêmes effectifs est aisée lorsque toutes les données présentent des coordonnées longitudinales x; différentes, une telle répartition peut s'avérer impossible lorsque plusieurs données partagent la même coordonnée longitudinale. La caractérisation des valeurs réelles et/ou simulées prises par la grandeur étudiée peut conduire à la même problématique. Notamment, lorsque les valeurs prises par la grandeur étudiée sont des valeurs discrètes, l'évaluation du modèle ne repose pas sur la comparaison des percentiles des valeurs simulées et observées, peu pertinents, mais sur la comparaison, entre les données simulées et les données observées, des probabilités d'obtenir les différentes valeurs discrètes possibles. Dans le cas de données pour lesquelles l'ensemble des valeurs possibles n'est pas borné, par exemple des données de comptage, il convient alors de regrouper les valeurs possibles en classes, puis d'estimer les probabilités de chaque classe. Cette méthode conduit à définir également des sous-intervalles de variation de la grandeur étudiée. Le choix de ces sous-intervalles est également crucial, puisque la caractérisation des valeurs réelles et/ou simulées en dépend directement. Le but de l'invention est donc de proposer un procédé d'analyse de données longitudinales ne présentant pas les inconvénients des procédés selon l'état de la technique, et permettant de caractériser de manière pertinente l'évolution de données. A cet effet, l'invention a pour objet un procédé d'analyse du type précité, caractérisé en ce que l'étape de détermination desdits sous-intervalles comprend : - la définition d'une fonction représentative d'une dispersion de ladite variable dans lesdits sous-intervalles, dont la valeur dépend des bornes inférieure et supérieure desdits sous-intervalles, et - la détermination des bornes inférieure et supérieure desdits sous-intervalles optimisant la valeur de ladite fonction. Selon d'autres aspects, le procédé d'analyse de données longitudinales comporte l'une ou plusieurs des caractéristiques suivantes : - ladite fonction dépend d'une somme des normes d'ordre p, avec p supérieur ou égal à 1, de la variable centrée sur lesdits sous-intervalles, - ladite fonction dépend d'une somme des variances de ladite variable sur lesdits sous-intervalles, - ladite fonction dépend en outre de la somme des variances des effectifs dans les différents sous-intervalles, - l'étape de détermination desdits sous-intervalles comprend la détermination de bornes inférieure et supérieure desdits sous-intervalles minimisant ladite fonction, - ladite fonction comprend un terme de pénalisation, croissant avec le nombre de sous-intervalles, - l'étape de détermination desdits sous-intervalles comprend en outre la détermination du nombre de sous-intervalles minimisant la valeur de ladite fonction,35 4 - ladite fonction comprend un terme pouvant s'exprimer sous la forme : K f = I l mi (zi - ak )p + f3Pen(KX ) k=1 i dans laquelle K, désigne le nombre de sous-intervalles, f3Pen(K,) est un terme de pénalisation, les termes zi désignent les valeurs prises par ladite variable sur le sous- intervalle d'indice k, et les termes mi désignent le nombre de répétitions de la valeur zi de ladite variable dans lesdites données. L'invention a également pour objet un programme d'ordinateur comportant des lignes de code qui lorsqu'elles sont exécutées par un calculateur, réalisent les étapes du procédé d'analyse selon l'invention, et un système d'analyse de données longitudinales, comprenant une unité de traitement apte à mettre en ceuvre le procédé selon l'invention, des moyens d'entrée de données longitudinales dans ladite unité de traitement, et une interface homme/machine comprenant des moyens d'affichage sous forme graphique desdites données. L'invention sera mieux comprise à l'aide de la description qui va suivre, donnée uniquement à titre d'exemple, et faite en se référant aux dessins annexés, sur lesquels : la figure 1 est une représentation graphique de données longitudinales, la figure 2 est un schéma illustrant un système d'analyse de données longitudinales selon un mode de réalisation de l'invention, la figure 3 est un schéma synoptique illustrant un procédé d'analyse de données longitudinales selon un mode de réalisation de l'invention, et - la figure 4 est une représentation graphique de données longitudinales telle qu'obtenue par le procédé d'analyse selon l'invention. La figure 1 est un exemple de représentation graphique de données longitudinales pouvant être analysées par le procédé selon l'invention. Ces données comprennent un ensemble de N couples (xi, yi ), i =1,2,...,N , xi et yi représentant les valeurs prises par deux variables aléatoires x et y, comprises dans deux intervalles notés respectivement lx et ly, bornés ou non bornés. Ces données représentent l'évolution de la variable y, c'est-à-dire de la grandeur étudiée, en fonction de la variable x, ou coordonnée longitudinale. On considérera par la suite que les coordonnées xi sont ordonnées par ordre croissant (i La détermination de sous-intervalles adjacents optimaux de variation de la variable x comprend ainsi la détermination du nombre KX de sous-intervalles et de l'ensemble ZT, c'est-à-dire des valeurs (2,,22,..., zK_,) des indices de la variable z définissant les bornes inférieures et supérieures de ces sous-intervalles. A cette fin, le procédé d'analyse comprend une étape 22 optionnelle de définition du nombre KX de sous-intervalles, lors de laquelle le nombre KX de sous-intervalles à déterminer est fixé arbitrairement par l'unité 12 de traitement ou par l'utilisateur. L'utilisateur peut cependant choisir de ne pas fixer le nombre KX de sous-intervalles à ce stade, auquel cas ce nombre KX sera déterminé de manière automatique et optimale dans la suite de l'analyse par l'unité 12 de traitement. La détermination de sous-intervalles adjacents optimisant la dispersion des données dans les différents sous-intervalles comprend ensuite une étape 24 de définition d'un critère d'optimalité du choix de ces sous-intervalles, c'est-à-dire d'une fonction JX représentative de cette dispersion, qui dépend du nombre KX de sous-intervalles et de leurs bornes inférieure et supérieure, c'est-à-dire du vecteur T. Cette fonction JX peut être définie par l'utilisateur ou choisie parmi plusieurs fonctions prédéfinies, en fonction des données étudiées. Elle s'exprime de manière générale sous la forme : JX (z, K) = F + f3Pen(K X ) dans laquelle F est une fonction caractérisant la dispersion des données dans les KX différents sous-intervalles, qui dépend des données contenues dans chacun des sous- intervalles, donc des bornes de ces intervalles, et /3Pen(KX) est une fonction croissante du nombre KX de sous-intervalles, appelée terme de pénalisation, /3 étant un paramètre pouvant être choisi par l'utilisateur. Kx La fonction F peut s'écrire sous la forme F = Fk , c'est-à-dire comme une k=1 somme de fonctions Fk mesurant chacune la dispersion des données dans le sous- intervalle d'indice k. La fonction Fk est par exemple définie comme une norme d'ordre p, avec p supérieur ou égal à 1, de la variable centrée sur le sous-intervalle Ik. Elle mesure alors comment les données sont réparties à l'intérieur du sous-intervalle Ik, et s'exprime sous la forme : ~k Fk = mi (zi - ak )p Ici, ak minimise Fk. Ainsi, si p=2, par exemple, alors ak est la moyenne pondérée Zk des valeurs prises par la variable z dans le sous-intervalle Ik, c'est-à-dire la moyenne pondérée des valeurs prises par la variable x dans le sous-intervalle Ik, définie par : La fonction F, notée F(», est ainsi égale à la somme des normes d'ordre p, avec p supérieur ou égal à 1, de la variable centrée sur les Kx sous-intervalles : Kx Kx rk F(1)-1Fk-~ lmi(zi - ak)p k=1 k=1 i=zk_1+1 Lorsque p est choisi égal à 2, la fonction Fk est alors proportionnelle à la variance des valeurs prises par la variable x sur le sous-intervalle Ik : rk Fk = 1 mi (zi - ak )2 La fonction F(» est alors égale, à un terme multiplicatif près, à la variance intra-15 intervalle de la variable x, Le. à la moyenne pondérée des variances de la variable x dans chaque sous-intervalle. La fonction F peut également être définie comme une fonction mesurant l'écart entre les effectifs des différents sous-intervalles Ik, c'est-à-dire le nombre de données pour lesquelles la variable x est comprise dans ce sous-intervalle Ik, et l'effectif moyen des 20 Kx sous-intervalles. La fonction F, notée Fe , s'exprime alors par : Kx / N F(2) =If nk-- Kx k=1 ~k où f désigne une fonction quelconque croissante, nk = 1 mi désigne l'effectif du sous- i=2k_ +1 intervalle d'indice k, et K représente l'effectif moyen des Ksous-intervalles. Kx La fonction F(2) est par exemple proportionnelle à la variance intra-intervalle des effectifs, et s'exprime par : Kx / 2 F(2) N = nk -- k=1 Kx La fonction F peut également être définie comme une combinaison linéaire des fonctions F(» et F(2), et mesure alors non seulement comment les données sont distribuées à l'intérieur de chaque sous-intervalle, mais aussi comment ces données sont distribuées entre les différents sous-intervalles. Le terme de pénalisation /3Pen(Kx) est choisi indépendamment de la fonction F. Il est par exemple proportionnel au nombre Kx de sous-intervalles : /3Pen(Kx) = /3Kx Le paramètre /3 peut dépendre du nombre N de données. Il est par exemple déterminé conformément à différentes approches de sélection de modèle, par exemple par minimisation d'un critère d'information tel que le critère d'information d'Aikake (AIC). Cependant, si le nombre Kx de sous-intervalles a été fixé lors de l'étape 22, /3 est choisi égal à zéro, de telle sorte que le terme de pénalisation /3Pen(Kx) est lui-même nul. La définition 24 de la fonction JX caractérisant la dispersion des données dans les différents sous-intervalles est suivie d'une étape 26 de détermination de bornes inférieure et supérieure des sous-intervalles optimisant la valeur de cette fonction JX. Lors de cette étape 26, l'unité 12 de traitement détermine ainsi le nombre Kx de sous-intervalles, s'il n'a pas été fixé lors de l'étape 22, et le vecteur 1- =(zz-1,zz-2,...,zz-K_I) minimisant la valeur de la fonctionJx(z,Kx). Cette étape peut être réalisée au moyen de tout type d'algorithmes de minimisation, par exemple suivant les algorithmes de programmation dynamique décrits dans le document « Using penalized contrasts for the change-point problem » (Lavielle M., Signal Processing, vol. 85, n. 8, pp 1501-1510, 2005). Si le nombre Kx de sous-intervalles a été fixé lors de l'étape 22, l'unité 12 de traitement détermine uniquement le vecteur z = (21,22,...,2K,H) minimisant la valeur de la fonction JX, c'est-à-dire la position des sous-intervalles optimisant la dispersion des données dans ces sous-intervalles. Notamment, si la fonction JX comprend un terme s'exprimant sous la forme F(» définie ci-dessus, c'est-à-dire dépend d'une somme des normes d'ordre p de la variable x centrée sur les sous-intervalles, le vecteur i déterminé lors de l'étape 26 est le vecteur optimisant la répartition des données au sein de chacun des sous-intervalles. Si la fonction JX comprend un terme s'exprimant sous la forme F(2) définie ci-dessus, le vecteur i optimise la distribution des données entre les différents sous-intervalles. Si le nombre KX de sous-intervalles n'a pas été fixé lors de l'étape 22, l'unité 12 de traitement détermine, outre le vecteurz = le nombre K, de sous-intervalles minimisant la valeur de la fonction JX, en établissant ainsi un compromis entre un grand nombre de sous-intervalles, souhaitable pour évaluer la variation de la grandeur étudiée, et un grand nombre de données par sous-intervalle, permettant de caractériser de manière plus précise les données au sein de chaque sous-intervalle. Ainsi, à l'issue de l'étape 26, les N données sont réparties en K, sous-intervalles adjacents, en fonction de leur coordonnée x;. Par exemple, une donnée de coordonnée x; telle que zzk Les Ky classes définissent Ky sous-intervalles adjacents de variation de la variable y. La détermination du nombre Ky de classes et de leurs bornes inférieure et supérieure est avantageusement réalisée de manière similaire à la détermination des K, sous-intervalles de variation de la variable x. Cette détermination comprend alors la définition d'une fonction Jy représentative de la dispersion de la variable y dans les Ky classes, dont la valeur dépend des bornes inférieure et supérieure de ces classes, et comprenant éventuellement un terme de pénalisation, et la détermination par l'unité 12 de traitement des bornes inférieure et supérieure de ces classes optimisant la valeur de la fonction Jy. Ce découpage automatique de l'intervalle de variation de la variable y en Ky classes optimise ainsi la distribution des N données dans les Ky classes. Ainsi, à l'issue de l'étape 28, chacun des K, sous-intervalles de variation de la variable x est associé à une ou plusieurs grandeurs Yk caractérisant les valeurs de la variable y sur ces sous-intervalles. Puis, lors d'une étape 30, l'unité 12 de traitement commande l'affichage sur les moyens 16 d'affichage des données analysées sous la forme d'un graphique, portant en abscisse la variable x, et en ordonnée la variable y, et représentant en outre les K, sous-intervalles déterminés lors de l'étape 26, ainsi que les paramètres Yk caractérisant les valeurs de la variable y sur ces sous-intervalles. On a ainsi illustré sur la figure 4 une représentation graphique de données longitudinales telles qu'affichées par les moyens d'affichage lors de l'étape 30. Sur ce graphique sont représentées les N données, sous forme de N points P;, de manière identique à la représentation illustrée en figure 1. Sont également représentés les K, sous-intervalles Ik de variation de la variable x, délimités par KX 1 lignes verticales Lk, obtenus en utilisant une fonction F(') telle que définie ci-dessus. Ce graphique représente en outre, sous forme de croix, les paramètres Yk caractérisant les valeurs de la variable y sur chaque sous-intervalle, qui sont ici les l e' 50à' et 90à' percentiles de la variable y, avec pour abscisse le paramètre ak défini dans la fonction F('). Ces percentiles sont reliés par des segments, de manière à visualiser leur évolution entre deux sous-intervalles Ik consécutifs. Le procédé d'analyse selon l'invention permet ainsi de déterminer de manière automatique, sans que l'intervention d'un expert soit nécessaire, des sous-intervalles optimaux de variation des variables étudiées, et donc d'avoir une évaluation qualitative et quantitative plus précise de la précision de modèles simulant des phénomènes réels. Cette optimalité revêt plusieurs aspects. Notamment, le procédé d'analyse permet d'optimiser à la fois la distribution des données à l'intérieur de chaque sous-intervalle, et la distribution des données entre les différents sous-intervalles, l'utilisateur restant libre de pondérer l'importance de ces deux critères dans la détermination des sous-intervalles. En outre, le procédé selon l'invention permet une détermination automatique et optimale du nombre de sous-intervalles, en établissant un compromis entre un nombre élevé de sous- intervalles favorisant une homogénéité des données au sein de chaque sous-intervalle et permettant de décrire plus précisément l'évolution de la grandeur étudiée, et un nombre élevé de données dans chaque sous-intervalle permettant une caractérisation plus précise de la grandeur étudiée. Il devra toutefois être compris que l'exemple de réalisation présenté ci-dessus n'est pas limitatif. Notamment, les sous-intervalles Ik déterminés lors de l'analyse des données observées peuvent être utilisés pour caractériser des données issues de simulations sur ces mêmes sous-intervalles. Les résultats de cette analyse sont alors avantageusement représentés sur le graphique de la figure 4, superposés aux résultats d'analyse des données observées. Une telle représentation permet ainsi à l'utilisateur de comparer les données observées aux données simulées, donc d'évaluer le modèle utilisé pour la simulation, par comparaison des paramètres Yk caractérisant les données observées et les données simulées sur chaque sous-intervalle. Par ailleurs, bien que le procédé selon l'invention ait été décrit précédemment dans le cadre de l'analyse de données de type VPC, il peut être appliqué à tout type de données longitudinales, caractérisant l'évolution d'au moins une première variable en fonction d'au moins une deuxième variable. REVENDICATIONS1.- Procédé d'analyse de données longitudinales caractérisant l'évolution d'au moins une première variable (y) relative à un phénomène réel en fonction d'au moins une deuxième variable (x), pour évaluer la précision de modèles simulant des phénomènes réels et comprenant des étapes de détermination (22, 24, 26) de sous-intervalles adjacents (Ik,lk) de variation d'au moins une desdites première et/ou deuxième variables (x, y) et de caractérisation (28) desdites données sur lesdits sous-intervalles (Ik, Ik ), caractérisé en ce que l'étape de détermination desdits sous-intervalles comprend : - la définition (24) d'une fonction (Jx, Jy) représentative d'une dispersion de ladite variable (x, y) dans lesdits sous-intervalles (Ik, lk ), dont la valeur dépend des bornes inférieure ) et supérieure (.zzk ) desdits sous-intervalles k,11,* ), et - la détermination (26) de bornes inférieure (z,k_i ) et supérieure (zJdesdits sous-intervalles (Ik,lk) optimisant la valeur de ladite fonction g, Jy). 2.- Procédé d'analyse selon la revendication 1, caractérisé en ce que ladite fonction (Jx, Jy) dépend d'une somme (FI des normes d'ordre p, avec p supérieur ou égal à 1, de la variable centrée sur lesdits sous-intervalles (Ik, Ik ). 3.- Procédé d'analyse selon l'une quelconque des revendications 1 ou 2, caractérisé en ce que ladite fonction (Jx, Jy) dépend d'une somme des variances de ladite variable sur lesdits sous-intervalles (Ik, lk ). 4.- Procédé d'analyse selon l'une quelconque des revendications précédentes, caractérisé en ce que ladite fonction (Jx, Jy) dépend en outre de la somme (F(2» des variances des effectifs (nx) dans les différents sous-intervalles (Ik, Ik ). 5.- Procédé d'analyse selon l'une quelconque des revendications précédentes, caractérisé en ce que l'étape (26) de détermination desdits sous-intervalles comprend la détermination de bornes inférieure et supérieure desdits sous-intervalles(,lk) minimisant ladite fonction (Jx, Jy). 6.- Procédé d'analyse selon l'une quelconque des revendications précédentes, caractérisé en ce que ladite fonction (Jx, Jy) comprend un terme de pénalisation (,6Pen(Kx) ), croissant avec le nombre (Kx) de sous-.intervalles. 7.- Procédé d'analyse selon la revendication 6, caractérisé en ce que l'étape (26) de détermination desdits sous-intervalle4,1k)comprend en outre la détermination du nombre (Kx) de sous-intervalles minimisant la valeur de ladite fonction (Jx, Jy). 8.- Procédé d'analyse selon l'une quelconque des revendications précédentes, caractérisé en ce que ladite fonction (JX, Jy) comprend un terme pouvant s'exprimer sous la forme : f = E E - ak )p + QPen(K, ) k=1 i dans laquelle K, désigne le nombre de sous-intervalles, flPen(Kjest un terme de pénalisation, les termes z; désignent les valeurs prises par ladite variable sur le sous-intervalle (Ik) d'indice k, et les termes m; désignent le nombre de répétitions de la valeur z; de ladite variable dans lesdites données. 9.- Programme d'ordinateur comportant des lignes de code qui lorsqu'elles sont exécutées par un calculateur, réalisent les étapes du procédé d'analyse selon l'une des revendications précédentes. 10.- Système (10) d'analyse de données longitudinales pour évaluer la précision de modèles simulant des phénomènes réels et comprenant une unité (12) de traitement apte à mettre en oeuvre le procédé selon l'une quelconque des revendications 1 à 8, des moyens (14) d'entrée de données longitudinales dans ladite unité (12) de traitement, et une interface homme/machine comprenant des moyens (16) d'affichage sous forme graphique desdites données.