Notes de cours Python scientifique

Numpy avancé

Pauli Virtanen

Traducteur

Jacques Beauregard

Numpy est à l'origine d'un ensemble d'outils scientifiques pour Python. Son but est simple : Mise en œuvre d'opérations efficaces sur de nombreux objets dans un bloc de mémoire. Comprendre comment cela fonctionne en détail permet d'en avoir une utilisation optimale : gestions de sa flexibilité, raccourcis utiles… et de construire un nouveau travail grâce à tout cela.

Ce tutoriel a pour but de couvrir :

- la structure des tableaux Numpy, et ce qui en découle. Trucs et astuces ;

- les fonctions universelles : quoi, pourquoi, et comment faire pour en créer une ;

- l'intégration avec d'autres outils: Numpy offre plusieurs façons d'encapsuler n'importe quelle donnée dans un ndarray, sans copies inutiles ;

- les fonctionnalités récemment ajoutées, et ce qu'il y a dedans selon moi: PEP 3118 tampons mémoires, fonctions universelles généralisées.

Retrouvez Les notes de cours scientifique de Python

Commentez

- Titre : Notes de cours Python scientifique

- Auteur : Pauli Virtanen

- Traducteur : Jacques Beauregard

- Parution : 26 avril 2018

- Licence : Les sources présentées sur cette page sont libres de droits et vous pouvez les utiliser à votre convenance. Par contre, la page de présentation constitue une œuvre intellectuelle protégée par les droits d'auteur. Copyright ® 2018 Pauli Virtanen. Aucune reproduction, même partielle, ne peut être faite de ce site et de l'ensemble de son contenu : textes, documents, images, etc. sans l'autorisation expresse de l'auteur. Sinon vous encourez selon la loi jusqu'à trois ans de prison et jusqu'à 300 000 € de dommages et intérêts.

Voir également :

- Notes de cours Python scientifique (Python Scientific Lecture Notes)

Sommaire

- I. Prérequis

- II. La vie d'un ndarray

- II-A. C'est…

- II-B. Bloc de mémoire

- II-C. Types de données

- II-C-1. Le descripteur

- II-C-2. Exemple : lire les fichiers .wav

- II-C-3. Casting et réinterprétation/vues

- II-C-4. Conversion de types

- II-C-4-a. Réinterprétation / Observation

- II-D. Plan d'indexation : les strides

- II-D-1. Point principal

- II-D-1-a. Programmation en C et Fortran

- II-D-1-b. Partage avec des nombres entiers

- II-D-2. Exemple : des dimensions faussées avec les strides

- II-D-3. Diffusion (Broadcasting)

- II-D-4. Plus de choses : diagonales

- II-D-5. Les effets de cache du CPU

- II-D-6. Exemple: les opérations en place (caveat emptor)

- II-E. Conclusions en bref

- III. Les fonctions universelles

- III-A. De quoi s'agit-il ?

- III-A-1. Les parties d'une fonction Ufunc

- III-A-2. Faciliter les choses

- III-B. Exercice : construire une ufunc à partir de zéro

- III-C. Solution : construire une ufunc à partir de zéro

- III-D. ufuncs généralisées

- IV. Les caractéristiques d'interopérabilité

- IV-A. Partager des données multidimensionnelles, typées

- IV-B. Le vieux protocole tampon

- IV-C. L'ancien protocole tampon

- IV-D. Protocole d'interface réseau

- V. Les tableaux fraternels : chararray, maskedarray, matrice

- V-A.  Chararray : opérations de chaînes de caractères vectorisées

- V-B. Données manquantes masked_array

- V-B-1. Le masque

- V-B-2. Fonctions du domaine-courant

- V-C.  Recarray : par pure convenance

- V-D. Matrix : convenance ?

- VI. En résumé

- VII. Contribuer à Numpy/Scipy

- VII-A. pourquoi

- VII-B. rapporter les erreurs

- VII-B-1. Exemple d'un bon rapport de bogue

I - Prérequis

- Numpy (>= 1.2; de préférence plus récent…)

- Cython (>= 0.12, pour l'exemple de la Fonction Universelle)

- PIL (utilisé dans deux ou trois exemples)

Dans cette partie, Numpy sera importé comme suit :

>>> import numpy as np

II - La vie d'un ndarray

II-A - C'est…

ndarray = bloc mémoire + schéma d'indexation + descripteur de type de données

- donnée brute

- comment localiser un élément

- comment interpréter un élément

1. typedef struct PyArrayObject {

2. PyObject_HEAD

4. /* Bloc de mémoire */

5. char *data;

7. /* Descripteur de type de données */

8. PyArray_Descr *descr;

10. /* Plan d'indexation */

11. int nd;

12. npy_intp *dimensions;

13. npy_intp *strides;

15. /* Autre chose */

16. PyObject *base;

17. int flags;

18. PyObject *weakreflist;

19. } PyArrayObject;

II-B - Bloc de mémoire

1. >>> x = np.array([1, 2, 3, 4], dtype=np.int32)

2. >>> x.data

3. <...at ...>

4. >>> str(x.data)

5. '\x01\x00\x00\x00\x02\x00\x00\x00\x03\x00\x00\x00\x04\x00\x00\x00'

Adresse mémoire de la donnée:

>>> x.__array_interface__['data'][0]
64803824

L'interface complète  _array_interface_ :

1. >>> x.__array_interface__

2. {'data': (35828928, False),

3. 'descr': [('', '<i4')],

4. 'shape': (4,),

5. 'strides': None,

6. 'typestr': '<i4',

7. 'version': 3}

deux ndarrays peuvent partager la même adresse mémoire

1. >>> x = np.array([1, 2, 3, 4])

2. >>> y = x[:-1]

3. >>> x[0] = 9

4. >>> y

5. array([9, 2, 3])

L'adresse mémoire n'a pas besoin d'appartenir à un ndarray :

>>> x = b'1234'
#Le 'b' est pour « bytes », nécessaire en Python 3

x est une chaîne de caractères (en Python 3, un « bytes »), nous pouvons le représenter comme un tableau d'entiers :

1. >>> y = np.frombuffer(x, dtype=np.int8)

2. >>> y.data

3. <... at ...>

4. >>> y.base is x

5. True

6. >>> y.flags

7. C_CONTIGUOUS : True

8. F_CONTIGUOUS : True

9. OWNDATA : False

10. WRITEABLE : False

11. ALIGNED : True

12. UPDATEIFCOPY : False

Les indicateurs owndata et writeable précisent le statut du bloc mémoire.

II-C - Types de données

II-C-1 - Le descripteur

dtype décrit un simple objet dans un tableau :

type

Type scalaire des données, parmi:
int8, int16, float64, et al. (taille fixée)
str, unicode, void (taille flexible)

itemsize

Taille du bloc de données

byteorder

Indicateur d'ordre des octets: big-endian > / little-endian < / non applicable |

fields

sub-dtypes, si c'est un type de données structurées

shape

Forme du tableau, si c'est un sous-tableau

1. >>> np.dtype(int).type

2. <type 'numpy.int64'>

3. >>> np.dtype(int).itemsize

4. 8

5. >>> np.dtype(int).byteorder

6. '='

II-C-2 - Exemple : lire les fichiers .wav

L'entête du fichier .wav :

chunk_id

"RIFF"

chunk_size

Nombre entier non signé à 4 octets au format little-endian

format

"WAVE"

fmt_id

"fmt "

fmt_size

Nombre entier non signé à 4 octets au format little-endian

audio_fmt

Nombre entier non signé à 2 octets au format little-endian

num_channels

Nombre entier non signé à 2 octets au format little-endian

sample_rate

Nombre entier non signé à 4 octets au format little-endian

byte_rate

Nombre entier non signé à 4 octets au format little-endian

block_align

Nombre entier non signé à 2 octets au format little-endian

bits_per_sample

Nombre entier non signé à 2 octets au format little-endian

data_id

"data"

data_size

Nombre entier non signé à 4 octets au format little-endian

- bloc de données brutes à 44 octets (au début du fichier)

- … suivi par  data_size d'octets de données réelles audio

L'entête du fichier .wav en tant que type de données structurées Numpy :

>>> wav_header_dtype = np.dtype([
...     ("chunk_id", (bytes, 4)),   # type scalaire à taille flexible, taille de l'objet 4
...     ("chunk_size", "<u4"),    # Nombre entier non signé à 32 bits au format little-endian
...     ("format", "S4"),         # chaîne de caractères à 4 octets
...     ("fmt_id", "S4"),
...     ("fmt_size", "<u4"),
...     ("audio_fmt", "<u2"),     #
...     ("num_channels", "<u2"),  # .. plus encore...
...     ("sample_rate", "<u4"),   #
...     ("byte_rate", "<u4"),
...     ("block_align", "<u2"),
...     ("bits_per_sample", "<u2"),
...     ("data_id", ("S1", (2, 2))), # sous-tableau, juste pour rire !
...     ("data_size", "u4"),
...     #
...     # Les données du son lui-même ne peuvent pas être représentées ici:
...     # ça n'a pas de taille fixe
...    ])

1. >>> wav_header_dtype['format']

2. dtype('|S4')

3. >>> wav_header_dtype.fields

4. <dictproxy object at ...>

5. >>> wav_header_dtype.fields['format']

6. (dtype('|S4'), 8)

- Le premier élément est le « sub-dtype » dans les données structurées, correspondant au nom format

- Le deuxième élément est son décalage ( en octets) à partir du début de l'objet

Exercice

Faire un dtype « creux » en utilisant les décalages, et seulement quelques champs :

1. >>> wav_header_dtype = np.dtype(dict(

2. ...   names=['format', 'sample_rate', 'data_id'],

3. ...   offsets=[offset_1, offset_2, offset_3], # compté à partir du début de la structure en octets

4. ...   formats=list of dtypes for each of the fields,

5. ... ))

et l'utiliser pour lire la fréquence d'échantillonnage et data_id (en tant que sous-tableau).

1. >>> f = open('data/test.wav', 'r')

2. >>> wav_header = np.fromfile(f, dtype=wav_header_dtype, count=1)

3. >>> f.close()

4. >>> print(wav_header)

5. [ ('RIFF', 17402L, 'WAVE', 'fmt ', 16L, 1, 1, 16000L, 32000L, 2, 16, [['d', 'a'], ['t', 'a']], 17366L)]

6. >>> wav_header['sample_rate']

7. array([16000], dtype=uint32)

Essayons d'accéder au sous-tableau :

1. >>> wav_header['data_id']

2. array([[['d', 'a'],

3. ['t', 'a']]],

4. dtype='|S1')

5. >>> wav_header.shape

6. (1,)

7. >>> wav_header['data_id'].shape

8. (1, 2, 2)

Quand on accède aux sous-tableaux, les dimensions sont ajoutées à la fin!

Il existe des modules comme wavfile, audiolab, etc. pour charger des données audio.

II-C-3 - Casting et réinterprétation/vues

Conversion de type

- sur attribution

- sur la construction d'un tableau

- sur l'arithmétique

- etc.

- et manuellement: .astype(dtype)

données réinterprétées

- manuellement: .view(dtype)

II-C-4 - Conversion de types

- En résumé, voici ce qu'il faut retenir du cast en arithmétique:

- seul le type (not value!) des opérandes nous intéresse

- le plus grand type « certifié » en mesure de représenter l'ensemble est considéré

- Les scalaires peuvent « perdre » face aux tableaux dans certaines situations

- Le cast copie les données en général:

1. >>> x = np.array([1, 2, 3, 4], dtype=np.float)

2. >>> x

3. array([ 1.,  2.,  3.,  4.])

4. >>> y = x.astype(np.int8)

5. >>> y

6. array([1, 2, 3, 4], dtype=int8)

7. >>> y + 1

8. array([2, 3, 4, 5], dtype=int8)

9. >>> y + 256

10. array([257, 258, 259, 260], dtype=int16)

11. >>> y + 256.0

12. array([ 257.,  258.,  259.,  260.])

13. >>> y + np.array([256], dtype=np.int32)

14. array([257, 258, 259, 260], dtype=int32)

- Cast sur un objet setté : le « dtype » du tableau n'est pas modifié sur l'affectation d'un objet :

1. >>> y[:] = y + 1.5

2. >>> y

3. array([2, 3, 4, 5], dtype=int8)

II-C-4-a - Réinterprétation / Observation

- Bloc de données en mémoire (4 octets)

0x01

||

0x02

||

0x03

||

0x04

- 4 entiers non signés uint8, ou,

- 4 entiers signés int8, ou,

- 2 entiers signés int16, ou,

- 1 entier signé int32, ou,

- 1 flottant float32, ou,

- …

- Comment passer de l'un à l'autre?

1. Changer le « dtype »:

1. >>> x = np.array([1, 2, 3, 4], dtype=np.uint8)

2. >>> x.dtype = "<i2"

3. >>> x

4. array([ 513, 1027], dtype=int16)

5. >>> 0x0201, 0x0403

6. (513, 1027)

0x01

0x02

||

0x03

0x04

Little-endian : l'octet le moins significatif est à gauche dans la mémoire

1. Créer une nouvelle vue :

1. >>> y = x.view("<i4")

2. >>> y

3. array([67305985], dtype=int32)

4. >>> 0x04030201

5. 67305985

0x01

0x02

0x03

0x04

Note

- .view crée des vues, il ne copie pas (ou n'altère pas) le bloc mémoire

- il ne change que le « dtype » (et ajuste la forme du tableau ):

1. >>> x[1] = 5

2. >>> y

3. array([328193], dtype=int32)

4. >>> y.base is x

5. True

Mini-exercice : réinterprétation des données

Vous avez des données RGBA dans un tableau:

>>> x = np.zeros((10, 10, 4), dtype=np.int8)
>>> x[:, :, 0] = 1
>>> x[:, :, 1] = 2
>>> x[:, :, 2] = 3
>>> x[:, :, 3] = 4

où les trois dernières dimensions sont les canaux R, B, et G, puis alpha.

Comment faire un tableau structuré (10, 10) avec des champs nommés ‘r', ‘g', ‘b', ‘a' sans copier les données ?

>>> y = ...

>>> assert (y['r'] == 1).all()
>>> assert (y['g'] == 2).all()
>>> assert (y['b'] == 3).all()
>>> assert (y['a'] == 4).all()

Solution

>>> y = x.view([('r', 'i1'),
... ('g', 'i1'),
... ('b', 'i1'),
... ('a', 'i1')]
... )[:, :, 0]

Un autre tableau qui prend exactement 4 octets de mémoire :

1. >>> y = np.array([[1, 3], [2, 4]], dtype=np.uint8).transpose()

2. >>> x = y.copy()

3. >>> x

4. array([[1, 2],

5. [3, 4]], dtype=uint8)

6. >>> y

7. array([[1, 2],

8. [3, 4]], dtype=uint8)

9. >>> x.view(np.int16)

10. array([[ 513],

11. [1027]], dtype=int16)

12. >>> 0x0201, 0x0403

13. (513, 1027)

14. >>> y.view(np.int16)

15. array([[ 769, 1026]], dtype=int16)

- Que s'est-il passé ?

- … nous avons besoin de regarder ce que x[0,1] vaut réellement

>>> 0x0301, 0x0402
(769, 1026)

II-D - Plan d'indexation : les strides

II-D-1 - Point principal

La question :

1. >>> x = np.array([[1, 2, 3],

2. ...              [4, 5, 6],

3. ...              [7, 8, 9]], dtype=np.int8)

4. >>> str(x.data)

5. '\x01\x02\x03\x04\x05\x06\x07\x08\t'

À quel octet débute l'objet x[1,2] dans x.data ?

La réponse (dans Numpy) :

- strides: le nombre d'octets à franchir pour trouver le prochain élément

- 1 stride par dimension

1. >>> x.strides

2. (3, 1)

3. >>> byte_offset = 3*1 + 1*2   # to find x[1,2]

4. >>> x.flat[byte_offset]

5. 6

6. >>> x[1, 2]

7. 6

- Simple, flexible

II-D-1-a - Programmation en C et Fortran

1. >>> x = np.array([[1, 2, 3],

2. ...               [4, 5, 6]], dtype=np.int16, order='C')

3. >>> x.strides

4. (6, 2)

5. >>> str(x.data)

6. '\x01\x00\x02\x00\x03\x00\x04\x00\x05\x00\x06\x00'

- Nécessité d'augmenter de 6 octets pour trouver la ligne suivante

- Nécessité d'augmenter de 2 octets pour trouver la colonne suivante

1. >>> y = np.array(x, order='F')

2. >>> y.strides

3. (2, 4)

4. >>> str(y.data)

5. '\x01\x00\x04\x00\x07\x00\x02\x00\x05\x00\x08\x00\x03\x00\x06\x00'

- Nécessité d'augmenter 2 octets pour trouver la ligne suivante

- Nécessité d'augmenter 4 octets pour trouver la colonne suivante

- De la même façon avec des dimensions supérieures

- C : les dernières dimensions varient plus vite (= strides plus petits)

- F : les premières dimensions varient plus vite

Maintenant nous pouvons comprendre le comportement de .view() :

>>> y = np.array([[1, 3], [2, 4]], dtype=np.uint8).transpose()
>>> x = y.copy()

La transposition n'affecte pas la mise en mémoire des données, seulement les strides

1. >>> x.strides

2. (2, 1)

3. >>> y.strides

4. (1, 2)

6. >>> str(x.data)

7. '\x01\x02\x03\x04'

8. >>> str(y.data)

9. '\x01\x03\x02\x04'

- les résultats sont différents lorsque c'est interprété comme deux « int16 »

- .copy() crée de nouveaux tableaux comme en C (par défaut)

II-D-1-b - Partage avec des nombres entiers

- Tout peut être représenté en changeant seulement shape, stride, et en ajustant au besoin le pointeur data !

- Ne jamais faire de copie des données

1. >>> x = np.array([1, 2, 3, 4, 5, 6], dtype=np.int32)

2. >>> y = x[::-1]

3. >>> y

4. array([6, 5, 4, 3, 2, 1], dtype=int32)

5. >>> y.strides

6. (-4,)

7. >>> y = x[2:]

8. >>> y.__array_interface__['data'][0] - x.__array_interface__['data'][0]

9. 8

10. >>> x = np.zeros((10, 10, 10), dtype=np.float)

11. >>> x.strides

12. (800, 80, 8)

13. >>> x[::2,::3,::4].strides

14. (1600, 240, 32)

De même, les transpositions ne créent jamais de copies (elles permutent juste les strides)

1. >>> x = np.zeros((10, 10, 10), dtype=np.float)

2. >>> x.strides

3. (800, 80, 8)

4. >>> x.T.strides

5. (8, 80, 800)

Mais toutes les opérations de redimensionnement de tableau ne peuvent être effectuées seulement en manipulant les strides :

1. >>> a = np.arange(6, dtype=np.int8).reshape(3, 2)

2. >>> b = a.T

3. >>> b.strides

4. (1, 2)

Jusqu'ici tout va bien, toutefois :

1. >>> str(a.data)

2. '\x00\x01\x02\x03\x04\x05'

3. >>> b

4. array([[0, 2, 4],

5. [1, 3, 5]], dtype=int8)

6. >>> c = b.reshape(3*2)

7. >>> c

8. array([0, 2, 4, 1, 3, 5], dtype=int8)

Ici, il n'y a aucun moyen de représenter le tableau c étant donné qu'il y a un stride et le bloc de mémoire pour a. Par conséquent, l'opération reshape a besoin d'en faire une copie.

II-D-2 - Exemple : des dimensions faussées avec les strides

Manipulation d'un stride

1. >>> from numpy.lib.stride_tricks import as_strided

2. >>> help(as_strided)

3. as_strided(x, shape=None, strides=None)

4. Fare un ndarray à partir du tableau avec les paramètres shape et strides fournis

as_strided ne vérifie pas que vous restez à l'intérieur des limites du bloc mémoire.

1. >>> x = np.array([1, 2, 3, 4], dtype=np.int16)

2. >>> as_strided(x, strides=(2*2, ), shape=(2, ))

3. array([1, 3], dtype=int16)

4. >>> x[::2]

5. array([1, 3], dtype=int16)

Exercice

1. array([1, 2, 3, 4], dtype=np.int8)

2. -> array([[1, 2, 3, 4],

3. [1, 2, 3, 4],

4. [1, 2, 3, 4]], dtype=np.int8)

en utilisant seulement as_strided :

Hint: byte_offset = stride[0]*index[0] + stride[1]*index[1] + ...

Solution

Stride peut aussi être 0

>>> x = np.array([1, 2, 3, 4], dtype=np.int8)
>>> y = as_strided(x, strides=(0, 1), shape=(3, 4))
>>> y
array([[1, 2, 3, 4],
[1, 2, 3, 4],
[1, 2, 3, 4]], dtype=int8)
>>> y.base.base is x
True

II-D-3 - Diffusion (Broadcasting)

- Faisons quelque chose d'utile avec ça : produit externe de [1, 2, 3, 4] et [5, 6, 7]

1. >>> x = np.array([1, 2, 3, 4], dtype=np.int16)

2. >>> x2 = as_strided(x, strides=(0, 1*2), shape=(3, 4))

3. >>> x2

4. array([[1, 2, 3, 4],

5. [1, 2, 3, 4],

6. [1, 2, 3, 4]], dtype=int16)

1. >>> y = np.array([5, 6, 7], dtype=np.int16)

2. >>> y2 = as_strided(y, strides=(1*2, 0), shape=(3, 4))

3. >>> y2

4. array([[5, 5, 5, 5],

5. [6, 6, 6, 6],

6. [7, 7, 7, 7]], dtype=int16)

1. >>> x2 * y2

2. array([[ 5, 10, 15, 20],

3. [ 6, 12, 18, 24],

4. [ 7, 14, 21, 28]], dtype=int16)

… semble quelque peu familier…

1. >>> x = np.array([1, 2, 3, 4], dtype=np.int16)

2. >>> y = np.array([5, 6, 7], dtype=np.int16)

3. >>> x[np.newaxis,:] * y[:,np.newaxis]

4. array([[ 5, 10, 15, 20],

5. [ 6, 12, 18, 24],

6. [ 7, 14, 21, 28]], dtype=int16)

- En interne, le principe de diffusion de tableau est en effet implémenté en utilisant des 0-strides.

II-D-4 - Plus de choses : diagonales

Défi

- Extraire la diagonale principale de la matrice : (en supposant que la mémoire est ordonnée comme en C):

1. >>> x = np.array([[1, 2, 3],

2. ...               [4, 5, 6],

3. ...               [7, 8, 9]], dtype=np.int32)

5. >>> x_diag = as_strided(x, shape=(3,), strides=(???,))

- Extraire les éléments de la première super-diagonale [2, 6].

- Et les sous-diagonales?

(Astuce pour les deux derniers points : le premier découpage déplace le point où le « striding » a commencé.)

Solution

Sélectionnez les diagonales:

>>> x_diag = as_strided(x, shape=(3, ), strides=((3+1)*x.itemsize, ))
>>> x_diag
array([1, 5, 9], dtype=int32)

Tranchez d'abord, pour ajuster le pointeur de données:

>>> as_strided(x[0, 1:], shape=(2, ), strides=((3+1)*x.itemsize, ))
array([2, 6], dtype=int32)

>>> as_strided(x[1:, 0], shape=(2, ), strides=((3+1)*x.itemsize, ))
array([4, 8], dtype=int32)

Note : Utilisation de np.diag

>>> y = np.diag(x, k=1)
>>> y
array([2, 6], dtype=int32)

Cependant,

>>> y.flags.owndata
False

Remarque : ce comportement a changé : avant numpy 1.9, np.diag en fait une copie.

Défi

Calculer la trace du tenseur :

1. >>> x = np.arange(5*5*5*5).reshape(5,5,5,5)

2. >>> s = 0

3. >>> for i in xrange(5):

4. ...    for j in xrange(5):

5. ...       s += x[j,i,j,i]

en utilisant les strides et sum() sur le résultat.

>>> y = as_strided(x, shape=(5, 5), strides=(TODO, TODO))
>>> s2 = ...
>>> assert s == s2

Solution

>>> y = as_strided(x, shape=(5, 5), strides=((5*5*5 + 5)*x.itemsize,
... (5*5 + 1)*x.itemsize))
>>> s2 = y.sum()

II-D-5 - Les effets de cache du CPU

La mise en mémoire peut affecter les performances:

1. In [1]: x = np.zeros((20000,))

3. In [2]: y = np.zeros((20000*67,))[::67]

5. In [3]: x.shape, y.shape

6. ((20000,), (20000,))

8. In [4]: %timeit x.sum()

9. 100000 loops, best of 3: 0.180 ms per loop

11. In [5]: %timeit y.sum()

12. 100000 loops, best of 3: 2.34 ms per loop

14. In [6]: x.strides, y.strides

15. ((8,), (536,))

Les strides plus petits sont-ils plus rapides ?

- Le CPU extrait des données de la mémoire principale à son cache en blocs

- Si de nombreux éléments du tableau sont ajustés en un seul bloc (petit stride):

- moins de transferts nécessaires

- plus rapide

Voir aussi : numexpr est conçu pour atténuer les effets de cache sur les traitements de tableaux.

II-D-6 - Exemple: les opérations en place (caveat emptor)

- Parfois,

>>> a -= b

n'est pas pareil que :

>>> a -= b.copy()

1. >>> x = np.array([[1, 2], [3, 4]])

2. >>> x -= x.transpose()

3. >>> x

4. array([[ 0, -1],

5. [ 4,  0]])

1. >>> y = np.array([[1, 2], [3, 4]])

2. >>> y -= y.T.copy()

3. >>> y

4. array([[ 0, -1],

5. [ 1,  0]])

- x et x.transpose() partagent les données

- x -= x.transpose() modifie les données élément par élément

- parce que x and x.transpose() ont un striding différent, des données modifiées réapparaissent sur le RHS

II-E - Conclusions en bref

- bloc mémoire : peut être partagé, .base, .data

- descripteur de type de données : données structurées, sous-tableaux, l'ordre des octets, cast, vues, .astype(), .view()

- indexation strided : les strides, programmation C-F, découpage w/ entiers, as_strided, diffusion, astuces sur le stride, diag, cohérence du cache du CPU

III - Les fonctions universelles

III-A - De quoi s'agit-il ?

- Les fonctions Ufunc (Universal functions) effectuent une opération élémentaire sur tous les éléments d'un tableau.

- Exemples:

np.add, np.subtract, scipy.special.*, ...

- elles supportent automatiquement : la diffusion, le cast.

- L'auteur d'une ufunc n'a qu'à fournir l'opération par élément, Numpy s'occupe du reste.

- L'opération par élément nécessite d'être implémenté en C (ou par exemple, Cython)

III-A-1 - Les parties d'une fonction Ufunc

1. Fournies par l'utilisateur

1. void ufunc_loop(void **args, int *dimensions, int *steps, void *data)

2. {

3. /*

4. * int8 output = elementwise_function(int8 input_1, int8 input_2)

5. *

6. * Cette fonction doit calculer la ufunc pour plusieurs valeurs à la fois

7. * de la manière suivante.

8. */

9. char *input_1 = (char*)args[0];

10. char *input_2 = (char*)args[1];

11. char *output = (char*)args[2];

12. int i;

14. for (i = 0; i < dimensions[0]; ++i) {

15. *output = elementwise_function(*input_1, *input_2);

16. input_1 += steps[0];

17. input_2 += steps[1];

18. output += steps[2];

19. }

20. }

2. La partie Numpy intégrée

1. char types[3]

3. types[0] = NPY_BYTE   /* type du premier argument en entrée */

4. types[1] = NPY_BYTE   /* type du second argument en entrée */

5. types[2] = NPY_BYTE   /* type du troisième argument en entrée */

7. PyObject *python_ufunc = PyUFunc_FromFuncAndData(

8. ufunc_loop,

9. NULL,

10. types,

11. 1, /* ntypes */

12. 2, /* num_inputs */

13. 1, /* num_outputs */

14. identity_element,

15. name,

16. docstring,

17. unused)

- Une ufunc peut également prendre en charge plusieurs combinaisons différentes de type entrée-sortie.

III-A-2 - Faciliter les choses

ufunc_loop est de forme très générique, et Numpy en fournit des préfabriquées

PyUfunc_f_f

float elementwise_func(float input_1)

PyUfunc_ff_f

float elementwise_func(float input_1, float input_2)

PyUfunc_d_d

double elementwise_func(double input_1)

PyUfunc_dd_d

double elementwise_func(double input_1, double input_2)

PyUfunc_D_D

elementwise_func(npy_cdouble *input, npy_cdouble* output)

PyUfunc_DD_D

elementwise_func(npy_cdouble *in1, npy_cdouble *in2, npy_cdouble* out)

- Seul elementwise_func doit être fourni

- … sauf si votre fonction par élément n'est pas dans l'une des formes ci-dessus.

III-B - Exercice : construire une ufunc à partir de zéro

La fractale de Mandelbrot est définie par l'itération où est un nombre complexe. Cette itération est répétée - si ne diverge pas, peu importe la durée du parcours des itérations, appartient à l'ensemble de Mandelbrot.

- Faire une ufunc appelée mandel(z0, c) qui calcule:

z = z0
for k in range(iterations):
z = z*z + c

Disons 100 itérations ou bien jusqu'à ce que z.real**2 + z.imag**2 > 1000. Utilisez-la pour déterminer quels éléments c sont dans l'ensemble de Mandelbrot.

- Notre fonction est simple, donc faire usage de l'aide de PyUFunc_*.

- L'écrire en Cython

Rappel : Quelques boucles ufunc préconstruites :

PyUfunc_f_f

float elementwise_func(float input_1)

PyUfunc_ff_f

float elementwise_func(float input_1, float input_2)

PyUfunc_d_d

double elementwise_func(double input_1)

PyUfunc_dd_d

double elementwise_func(double input_1, double input_2)

PyUfunc_D_D

elementwise_func(complex_double *input, complex_double* output)

PyUfunc_DD_D

elementwise_func(complex_double *in1, complex_double *in2, complex_double* out)

Codes types :

NPY_BOOL, NPY_BYTE, NPY_UBYTE, NPY_SHORT, NPY_USHORT, NPY_INT, NPY_UINT,
NPY_LONG, NPY_ULONG, NPY_LONGLONG, NPY_ULONGLONG, NPY_FLOAT, NPY_DOUBLE,
NPY_LONGDOUBLE, NPY_CFLOAT, NPY_CDOUBLE, NPY_CLONGDOUBLE, NPY_DATETIME,
NPY_TIMEDELTA, NPY_OBJECT, NPY_STRING, NPY_UNICODE, NPY_VOID

III-C - Solution : construire une ufunc à partir de zéro

1. # La fonction par élément

2. # ------------------------

4. cdef void mandel_single_point(double complex *z_in,

5. double complex *c_in,

6. double complex *z_out) nogil:

7. #

8. # L'itération de Mandelbrot

9. #

11. #

12. # Quelques points:

13. #

14. # - Il n'est pas autorisé d'appeler n'importe qu'elles

15. #   fonctions Python ici.

16. #

17. #   La boucle Ufunc fonctionne avec le Python Global Interpreter Lock publié

18. #   Par conséquent, le ``nogil``.

19. #

20. # - Et toutes les variables locales doivent être déclarées avec ''cdef''

21. #

22. #  Notez également que cette fonction reçoit des *pointeurs* vers les données

23. #   La solution « traditionnelle » pour passer des variables complexes autour

24. #

26. cdef double complex z = z_in[0]

27. cdef double complex c = c_in[0]

28. cdef int k  # L'entier utilisé dans la boucle for

30. # Itération simple

32. for k in range(100):

33. z = z*z + c

34. if z.real**2 + z.imag**2 > 1000:

35. break

37. # Renvoie la réponse pour ce point

38. z_out[0] = z

41. # Définitions Boilerplate Cython

42. #

43. # Vous n'avez pas nécessairement besoin de lire cette partie, c'est simplement extrait des entêtes de Numpy C.

44. # ----------------------------------------------------------

46. cdef extern from "numpy/arrayobject.h":

47. void import_array()

48. ctypedef int npy_intp

49. cdef enum NPY_TYPES:

50. NPY_CDOUBLE

52. cdef extern from "numpy/ufuncobject.h":

53. void import_ufunc()

54. ctypedef void (*PyUFuncGenericFunction)(char**, npy_intp*, npy_intp*, void*)

55. object PyUFunc_FromFuncAndData(PyUFuncGenericFunction* func, void** data,

56. char* types, int ntypes, int nin, int nout,

57. int identity, char* name, char* doc, int c)

59. void PyUFunc_DD_D(char**, npy_intp*, npy_intp*, void*)

62. # module initialisation requis

63. # ------------------------------

65. import_array()

66. import_ufunc()

69. # La déclaration ufunc réelle

70. # ----------------------------

72. cdef PyUFuncGenericFunction loop_func[1]

73. cdef char input_output_types[3]

74. cdef void *elementwise_funcs[1]

76. loop_func[0] = PyUFunc_DD_D

78. input_output_types[0] = NPY_CDOUBLE

79. input_output_types[1] = NPY_CDOUBLE

80. input_output_types[2] = NPY_CDOUBLE

82. elementwise_funcs[0] = <void*>mandel_single_point

84. mandel = PyUFunc_FromFuncAndData(

85. loop_func,

86. elementwise_funcs,

87. input_output_types,

88. 1, # nombre de type d'entrées supportées

89. 2, # nombre d'arguments en entrée

90. 1, # nombre d'arguments en sortie

91. 0, # élément `identity`, ne vous en occupez jamais

92. "mandel", # nom de fonction

93. "mandel(z, c) -> computes iterated z*z + c", # docstring

94. 0 # inutilisé

95. )

96. import numpy as np

97. import mandel

98. x = np.linspace(-1.7, 0.6, 1000)

99. y = np.linspace(-1.4, 1.4, 1000)

100. c = x[None,:] + 1j*y[:,None]

101. z = mandel.mandel(c, c)

103. import matplotlib.pyplot as plt

104. plt.imshow(abs(z)**2 < 1000, extent=[-1.7, 0.6, -1.4, 1.4])

105. plt.gray()

106. plt.show()

La plupart des outils pourraient être automatisés par ces modules Cython.

Plusieurs types d'entrées acceptés

Par exemple supportant les versions à simple et double précision

1. cdef void mandel_single_point(double complex *z_in,

2. double complex *c_in,

3. double complex *z_out) nogil:

4. ...

6. cdef void mandel_single_point_singleprec(float complex *z_in,

7. float complex *c_in,

8. float complex *z_out) nogil:

9. ...

11. cdef PyUFuncGenericFunction loop_funcs[2]

12. cdef char input_output_types[3*2]

13. cdef void *elementwise_funcs[1*2]

15. loop_funcs[0] = PyUFunc_DD_D

16. input_output_types[0] = NPY_CDOUBLE

17. input_output_types[1] = NPY_CDOUBLE

18. input_output_types[2] = NPY_CDOUBLE

19. elementwise_funcs[0] = <void*>mandel_single_point

21. loop_funcs[1] = PyUFunc_FF_F

22. input_output_types[3] = NPY_CFLOAT

23. input_output_types[4] = NPY_CFLOAT

24. input_output_types[5] = NPY_CFLOAT

25. elementwise_funcs[1] = <void*>mandel_single_point_singleprec

27. mandel = PyUFunc_FromFuncAndData(

28. loop_func,

29. elementwise_funcs,

30. input_output_types,

31. 2, # nombre de types d'entrées supportées   <-------------

32. 2, # nombre d'arguments en entrée

33. 1, # nombre d'arguments en sortie

34. 0, # élément `identity`, ne vous en occupez jamais

35. "mandel", # nom de fonction

36. "mandel(z, c) -> computes iterated z*z + c", # docstring

37. 0 # inutilisé

38. )

III-D - ufuncs généralisées

ufunc

output = elementwise_function(output)

output et input peuvent être tous les deux un élément d'un tableau

ufunc généralisée

output et input peuvent être des tableaux avec un nombre de dimensions fixe.
Par exemple, la matrice de trace (somme des éléments diagonaux)

1. input shape = (n, n)

2. output shape = ()      i.e.  scalar

4. (n, n) -> ()

Produit matriciel :

1. input_1 shape = (m, n)

2. input_2 shape = (n, p)

3. output shape  = (m, p)

5. (m, n), (n, p) -> (m, p)

- Ceci est appelé la « signature » de la fonction ufunc généralisée

- Les dimensions sur lesquelles le g-ufunc agit, sont des « dimensions fondamentales »

Statuts dans Numpy

- g-ufuncs sont déjà en Numpy …

- Les nouvelles peuvent être créées avec PyUFunc_FromFuncAndDataAndSignature

- … Mais nous ne livrons pas avec des g-ufuncs publiques, sauf pour les tests, ATM

1. >>> import numpy.core.umath_tests as ut

2. >>> ut.matrix_multiply.signature

3. '(m,n),(n,p)->(m,p)'

>>> x = np.ones((10, 2, 4))
>>> y = np.ones((10, 4, 5))
>>> ut.matrix_multiply(x, y).shape
(10, 2, 5)

- Les deux dernières dimensions sont devenues des dimensions fondamentales, et sont modifiées tout comme la signature

- Par ailleurs, la g-func fonctionne « par élément »

- La multiplication matricielle pourrait de cette façon être utile pour fonctionner sur de nombreuses petites matrices à la fois.

Boucle ufunc généralisée

Multiplication matricielle (m,n),(n,p) -> (m,p)

1. void gufunc_loop(void **args, int *dimensions, int *steps, void *data)

2. {

3. char *input_1 = (char*)args[0];  /* idem que précédemment */

4. char *input_2 = (char*)args[1];

5. char *output = (char*)args[2];

7. int input_1_stride_m = steps[3];  /* Les strides pour les dimensions fondamentales*/

8. int input_1_stride_n = steps[4];  /* sont ajoutés après la partie non-fondamentale */

9. int input_2_strides_n = steps[5]; /* étapes */

10. int input_2_strides_p = steps[6];

11. int output_strides_n = steps[7];

12. int output_strides_p = steps[8];

14. int m = dimension[1]; /* les dimensions fondamentales sont ajoutées après */

15. int n = dimension[2]; /* la dimension principale  */

16. int p = dimension[3]; /* signature */

18. int i;

20. for (i = 0; i < dimensions[0]; ++i) {

21. matmul_for_strided_matrices(input_1, input_2, output,

22. strides for each array...);

24. input_1 += steps[0];

25. input_2 += steps[1];

26. output += steps[2];

27. }

28. }

IV - Les caractéristiques d'interopérabilité

IV-A - Partager des données multidimensionnelles, typées

Supposons que vous :

1. Écriviez une bibliothèque qui gère des données binaires multidimensionnelles,

2. Vous vouliez faciliter la manipulation des données avec Numpy, ou une tout autre bibliothèque,

3. … mais que vous ne vouliez pas avoir Numpy comme dépendance.

Actuellement, trois solutions :

1. la “vieille” interface tampon

2. l'interface réseau

3. la “nouvelle” interface tampon (PEP 3118)

IV-B - Le vieux protocole tampon

- Uniquement des tampons 1-D

- Pas d'information concernant le type de données

- Interface C-level ; PyBufferProcs tp_as_buffer dans l'objet type

- Mais c'est intégré dans Python (par exemple les chaînes de caractères le supportent)

Miniexercice utilisant PIL (Python Imaging Library):

1. >>> from PIL import Image

2. >>> data = np.zeros((200, 200, 4), dtype=np.int8)

3. >>> data[:, :] = [255, 0, 0, 255] # Red

4. >>> # Dans PIL, les images RGBA consistent en des entiers de 32-bit dont les octets sont [RR,GG,BB,AA]

5. >>> data = data.view(np.int32).squeeze()

6. >>> img = Image.frombuffer("RGBA", (200, 200), data)

7. >>> img.save('test.png')

Q: Vérifiez si les données sont à présent modifiées, et que l'img soit encore sauvegardée.

IV-C - L'ancien protocole tampon

1. " "

2. De la  mémoire tampon

3. ============

5. Montrer comment échanger des données entre numpy et une bibliothèque qui ne connaît que l'interface tampon.

7. import numpy as np

8. import Image

10. # faisons une image simple, au format RGBA

12. x = np.zeros((200, 200, 4), dtype=np.int8)

14. x[:,:,0] = 254 # rouge

15. x[:,:,3] = 255 # opaque

17. data = x.view(np.int32) # Vérifiez que vous comprenez pourquoi cela est OK !

19. img = Image.frombuffer("RGBA", (200, 200), data)

20. img.save('test.png')

22. #

23. # Modifiez les données d'origine, et enregistrez à nouveau

24. #

25. # Il se trouve que PIL, qui ne connaît presque rien sur Numpy,

26. # partage heureusement les mêmes données

27. #

29. x[:,:,1] = 254

30. img.save('test2.png')

IV-D - Protocole d'interface réseau

- Tampons multidimensionnels

- Présence d'information pour chaque type de données

- Approche spécifique Numpy; lentement déprécié (mais ne va pas disparaître)

- Non intégré dans Python autrement

- Voir aussi : Documentation: http://docs.scipy.org/doc/numpy/reference/arrays.interface.html

1. >>> x = np.array([[1, 2], [3, 4]])

2. >>> x.__array_interface__

3. {'data': (171694552, False),      # adresse mémoire de la donnée, en lecture seule ?

4. 'descr': [('', '<i4')],          # descripteur du type de donnée

5. 'typestr': '<i4',                # idem, dans une autre forme

6. 'strides': None,                 # strides; ou rien si c'est du C

7. 'shape': (2, 2),

8. 'version': 3,

9. }

1. >>> import Image

2. >>> img = Image.open('data/test.png')

3. >>> img.__array_interface__

4. {'data': ...,

5. 'shape': (200, 200, 4),

6. 'typestr': '|u1'}

7. >>> x = np.asarray(img)

8. >>> x.shape

9. (200, 200, 4)

10. >>> x.dtype

11. dtype('uint8')

Note : Une variante proche du C de l'interface réseau est également définie

V - Les tableaux fraternels : chararray, maskedarray, matrice

V-A -  Chararray : opérations de chaînes de caractères vectorisées

1. >>> x = np.array(['a', '  bbb', '  ccc']).view(np.chararray)

2. >>> x.lstrip(' ')

3. chararray(['a', 'bbb', 'ccc'],

4. dtype='|S5')

5. >>> x.upper()

6. chararray(['A', '  BBB', '  CCC'],

7. dtype='|S5')

.view() a un autre sens : il peut faire d'un ndarray une instance d'une sous-classe d'un ndarray spécialisé

V-B - Données manquantes masked_array

Les tableaux masqués sont des tableaux qui peuvent avoir des entrées manquantes ou invalides. Par exemple, supposons que nous ayons un tableau où la quatrième entrée est invalide :

>>> x = np.array([1, 2, 3, -99, 5])

Une façon de décrire cela est de créer un tableau masqué

1. >>> mx = np.ma.masked_array(x, mask=[0, 0, 0, 1, 0])

2. >>> mx

3. masked_array(data = [1 2 3 -- 5],

4. mask = [False False False  True False],

5. fill_value = 999999)

Le mean masqué ignore les données masquées

>>> mx.mean()
2.75
>>> np.mean(mx)
2.75

Toutes les fonctions Numpy ne respectent pas les masques, par exemple np.dot, donc vérifier les types de retour

Le masked_array renvoie une vue au tableau d'origine :

>>> mx[1] = 9
>>> x
array([  1,   9,   3, -99,   5])

V-B-1 - Le masque

Vous pouvez modifier le masque par assignation :

1. >>> mx[1] = np.ma.masked

2. >>> mx

3. masked_array(data = [1 -- 3 -- 5],

4. mask = [False  True False  True False],

5. fill_value = 999999)

Le masque est ainsi effacé :

1. >>> mx[1] = 9

2. >>> mx

3. masked_array(data = [1 9 3 -- 5],

4. mask = [False False False  True False],

5. fill_value = 999999)

Le masque est également disponible directement :

>>> mx.mask
array([False, False, False,  True, False], dtype=bool)

Les entrées masquées peuvent être remplies avec une valeur déterminée pour obtenir un tableau en retour :

1. >>> x2 = mx.filled(-1)

2. >>> x2

3. array([ 1,  9,  3, -1,  5])

Le masque peut aussi être effacé

1. >>> mx.mask = np.ma.nomask

2. >>> mx

3. masked_array(data = [1 9 3 -99 5],

4. mask = [False False False False False],

5. fill_value = 999999)

V-B-2 - Fonctions du domaine-courant

Le package de tableau masqué contient également des fonctions de domaine-courant :

1. >>> np.ma.log(np.array([1, 2, -1, -2, 3, -5]))

2. masked_array(data = [0.0 0.69314718056 -- -- 1.09861228867 --],

3. mask = [False False  True  True False  True],

4. fill_value = 1e+20)

La Prise en charge simplifiée est plus transparente pour traiter les données manquantes dans les tableaux dans Numpy 1.7. Restez à l'écoute !

Exemple : Statistiques masquées

Des rangers canadiens ont été distraits lors du comptage de lièvres et de lynx en 1903-1910 et 1917-1918, et ont obtenu des chiffres incorrects. (Des cultivateurs de carottes sont restés en alerte, cependant.) Calculez les populations moyennes au fil du temps, en ignorant les chiffres invalides.

1. >>> data = np.loadtxt('data/populations.txt')

2. >>> populations = np.ma.masked_array(data[:,1:])

3. >>> year = data[:, 0]

5. >>> bad_years = (((year >= 1903) & (year <= 1910))

6. ...            | ((year >= 1917) & (year <= 1918)))

7. >>> # '&' signifie 'and' et '|' signifie 'or'

8. >>> populations[bad_years, 0] = np.ma.masked

9. >>> populations[bad_years, 1] = np.ma.masked

11. >>> populations.mean(axis=0)

12. masked_array(data = [40472.7272727 18627.2727273 42400.0],

13. mask = [False False False],

14. fill_value = 1e+20)

16. >>> populations.std(axis=0)

17. masked_array(data = [21087.656489 15625.7998142 3322.50622558],

18. mask = [False False False],

19. fill_value = 1e+20)

Notez que Matplotlib connaît les tableaux masqués.

>>> plt.plot(year, populations, 'o-')
[<matplotlib.lines.Line2D object at ...>, ...]

V-C -  Recarray : par pure convenance

1. >>> arr = np.array([('a', 1), ('b', 2)], dtype=[('x', 'S1'), ('y', int)])

2. >>> arr2 = arr.view(np.recarray)

3. >>> arr2.x

4. chararray(['a', 'b'],

5. dtype='|S1')

6. >>> arr2.y

7. array([1, 2])

V-D - Matrix : convenance ?

- Toujours en 2D

- * est le produit matriciel, et non la matrice par élément

>>> np.matrix([[1, 0], [0, 1]]) * np.matrix([[1, 2], [3, 4]])
matrix([[1, 2],
[3, 4]])

VI - En résumé

- Structure des ndarray : données, dtypes, strides.

- Fonctions universelles: opérations par élément, comment faire pour en créer une ?

- Sous-classes ndarray

- Différentes interfaces de tampons pour l'intégration avec d'autres outils

- Ajouts récents : PEP 3118, ufuncs généralisées

VII - Contribuer à Numpy/Scipy

Récupérez ce tutoriel: http://www.euroscipy.org/talk/882

VII-A - pourquoi

- Pour signaler un bogue

- Pour mieux comprendre

- Pour obtenir un code de qualité

- Pour apporter son aide

VII-B - rapporter les erreurs

- Bug tracker (préférez cela)

- http://www.numpy.org/

- https://scipy.org/scipylib/index.html

- Cliquez sur le lien “S'inscrire” pour obtenir un compte

- listes de diffusion ( scipy.org/Mailing_Lists )

- Si vous n'êtes pas sûr

- Aucune réponse au bout d'une semaine ou deux? Déposez juste un ticket d'erreur

VII-B-1 - Exemple d'un bon rapport de bogue

Titre : numpy.random.permutations échoue pour des arguments de type non entier

Je suis en train de générer des permutations aléatoires, en utilisant numpy.random.permutations

Lors de l'appel de numpy.random.permutation avec des arguments non entiers.

Cela échoue avec un message d'erreur sibyllin :

>>> np.random.permutation(12)
array([ 6, 11,  4, 10,  2,  8,  1,  7,  9,  3,  0,  5])
>>> np.random.permutation(12.)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "mtrand.pyx", line 3311, in mtrand.RandomState.permutation
File "mtrand.pyx", line 3254, in mtrand.RandomState.shuffle
TypeError: len() of unsized object

Cela se produit aussi avec des arguments de type long, et donc

np.random.permutation(X.shape[0]) où X est un tableau qui échoue sur Windows 64-bits(où la forme est un tuple de longs).

Ce serait génial si ça pouvait caster en entier ou au moins lever une erreur spécifique pour des types non entiers.

Je suis en train d'utiliser Numpy 1.4.1, construit à partir de l'archive officielle, sur Windows

64 avec Visual studio 2008, sur Python.org 64-bit Python.

1. Qu'êtes-vous en train de faire?

2. Extrait de code reproduisant l'erreur (si possible)

- Que se passe-t-il réellement ?

- Qu'est-ce que vous attendez ?

3. Plate-forme(Windows / Linux / OSX, 32/64 bits, x86/PPC…)

4. Version de Numpy/Scipy

>>> print np.__version__
2...

5. Vérifiez que ce qui suit est ce que vous aviez prévu :

>>> print np.__file__
/...

6. Dans le cas où vous traîniez de vieilles installations de Numpy ou endommagées.
En cas de doute, essayez de supprimer les installations existantes de Numpy, puis réinstallez…