[PYTHON] Analyse des données à l'aide de xarray

Analyse des données à l'aide de xarray

Article précédent

Je présenterai l'analyse des données à l'aide de xarray présentée ci-dessus d'un point de vue plus pratique.

A propos de l'auteur

Depuis 2017, il est membre du développement de pydata xarray. Récemment, je suis principalement un utilisateur, mais j'envoie occasionnellement des PR pour des corrections de bogues.

J'utilise xarray pour presque toutes les analyses de données que je fais dans mon entreprise principale.

À propos du lecteur prévu

On suppose que vous avez une certaine connaissance de numpy (vous pouvez effectuer des opérations et des opérations en utilisant np.ndarray).

Structure de cet article

Dans cet article, je présenterai à nouveau le modèle de données xarray, puis discuterai du fonctionnement de base de l'indexation. Il peut être prudent de dire que cette indexation est l'essence même de la radiographie.

Après cela, je présenterai concat, qui est l'opération inverse de l'indexation.

Sur cette base, lors de l'examen de la création réelle de données, je présenterai le type de structure de données à utiliser pour faciliter l'analyse avec xarray.

Je parlerai également un peu de la sauvegarde et du chargement des données et du traitement hors mémoire à l'aide de dask.

Cet article a été rédigé à l'aide de Google Collaboratory. Tout le monde peut exécuter et tracer en copiant depuis ici.

modèle de données xarray

Aperçu

En un mot, xarray est ** un tableau multidimensionnel avec des axes **. Vous pouvez le considérer comme une combinaison de nd-array de numpy et de pd.Series de pandas.

Une fois que vous vous y êtes habitué, c'est très pratique car vous n'avez pas à penser aux coordonnées lors de l'analyse des données.

Par contre, on dit souvent que l'utilisation est un peu unique et que le coût d'apprentissage est élevé. Cet article vise à réduire autant que possible ce coût d'apprentissage.

Il existe deux classes principales dans xarray. L'un est xarray.DataArray et l'autre est xarray.Dataset. Je vais les présenter.

C'est presque le même que le contenu décrit dans Multidimensional data analysis library xarray.

xarray.DataArray

Présentation de xarray.DataArray, un ** tableau multidimensionnel avec axes **. J'expliquerai comment créer un objet plus tard, mais j'expliquerai d'abord l'utilisation des données dans xarray.tutorial.

import numpy as np
import xarray as xr
data = xr.tutorial.open_dataset('air_temperature')['air']
data
<xarray.DataArray 'air' (time: 2920, lat: 25, lon: 53)>
[3869000 values with dtype=float32]
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

data est une instance 3D xr.DataArray. L'abréviation officielle de xarray est «xr».

Lorsque vous affichez data, vous verrez un résumé comme celui ci-dessus. Le (time: 2920, lat: 25, lon: 53) dans la rangée du haut représente le nom de dimension du tableau stocké et son nombre d'éléments.

Comme vous pouvez le voir, xarray ** traite chaque dimension avec un nom ** (appelé * label * dans la documentation officielle).

Le premier axe est l'axe «temps», le deuxième axe est l'axe «lat» et le troisième axe est l'axe «lon». En le nommant ainsi, vous n'avez pas à vous soucier de l'ordre des axes.

De plus, des valeurs de coordonnées sont ajoutées à chaque axe. Une liste de coordonnées est affichée dans la section «Coordonnées». Dans ces données, tous les axes de «temps», «lat», «lon» ont des coordonnées.

Vous pouvez également joindre d'autres données. Ils sont dans la section «Attributs».

Comme décrit en détail dans la partie indexation, les axes de coordonnées permettent d'obtenir intuitivement des données dans la plage souhaitée. Par exemple, si vous voulez les données du 30 août 2014,

data.sel(time='2014-08-30')
<xarray.DataArray 'air' (time: 4, lat: 25, lon: 53)>
array([[[270.29   , 270.29   , ..., 258.9    , 259.5    ],
        [273.     , 273.1    , ..., 262.6    , 264.19998],
        ...,
        [299.5    , 298.79   , ..., 299.     , 298.5    ],
        [299.29   , 299.29   , ..., 299.9    , 300.     ]],

       [[269.9    , 270.1    , ..., 258.6    , 259.5    ],
        [273.     , 273.29   , ..., 262.5    , 265.     ],
        ...,
        [299.19998, 298.9    , ..., 298.5    , 298.     ],
        [299.4    , 299.69998, ..., 299.19998, 299.6    ]],

       [[270.4    , 270.69998, ..., 261.29   , 261.9    ],
        [273.     , 273.6    , ..., 266.4    , 268.6    ],
        ...,
        [297.9    , 297.6    , ..., 298.29   , 298.29   ],
        [298.69998, 298.69998, ..., 299.1    , 299.4    ]],

       [[270.     , 270.4    , ..., 266.     , 266.19998],
        [273.     , 273.9    , ..., 268.1    , 269.69998],
        ...,
        [298.5    , 298.29   , ..., 298.69998, 299.     ],
        [299.1    , 299.19998, ..., 299.5    , 299.69998]]], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2014-08-30 ... 2014-08-30T18:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Vous pouvez le faire comme ça. Ici, la méthode .sel est utilisée pour sélectionner des données en utilisant les axes de coordonnées. .sel (time = '2014-08-30') équivaut à référencer l'axe temps et à sélectionner les données 2014-08-30 ''.

Opérations sur DataArray

DataArray est également un tableau multidimensionnel comme np.ndarray, vous pouvez donc effectuer des opérations comme sur np.ndarray.

data * 2
<xarray.DataArray 'air' (time: 2920, lat: 25, lon: 53)>
array([[[482.4    , 485.     , 487.     , ..., 465.59998, 471.     ,
         477.19998],
        [487.59998, 489.     , 489.4    , ..., 465.59998, 470.59998,
         478.59998],
        [500.     , 499.59998, 497.78   , ..., 466.4    , 472.78   ,
         483.4    ],
        ...,
        [593.2    , 592.39996, 592.8    , ..., 590.8    , 590.2    ,
         589.39996],
        [591.8    , 592.39996, 593.58   , ..., 591.8    , 591.8    ,
         590.39996],
        [592.58   , 593.58   , 594.2    , ..., 593.8    , 593.58   ,
         593.2    ]],

       [[484.19998, 485.4    , 486.19998, ..., 464.     , 467.19998,
         471.59998],
        [487.19998, 488.19998, 488.4    , ..., 462.     , 465.     ,
         471.4    ],
        [506.4    , 505.78   , 504.19998, ..., 461.59998, 466.78   ,
         477.     ],
...
        [587.38   , 587.77997, 590.77997, ..., 590.18   , 589.38   ,
         588.58   ],
        [592.58   , 594.38   , 595.18   , ..., 590.58   , 590.18   ,
         588.77997],
        [595.58   , 596.77997, 596.98   , ..., 591.38   , 590.98   ,
         590.38   ]],

       [[490.18   , 488.58   , 486.58   , ..., 483.37997, 482.97998,
         483.58   ],
        [499.78   , 498.58   , 496.78   , ..., 479.18   , 480.58   ,
         483.37997],
        [525.98   , 524.38   , 522.77997, ..., 479.78   , 485.18   ,
         492.58   ],
        ...,
        [587.58   , 587.38   , 590.18   , ..., 590.58   , 590.18   ,
         589.38   ],
        [592.18   , 593.77997, 594.38   , ..., 591.38   , 591.38   ,
         590.38   ],
        [595.38   , 596.18   , 596.18   , ..., 592.98   , 592.38   ,
         591.38   ]]], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00

Comme cela devrait être intuitivement, la multiplication est effectuée uniquement sur les données du tableau et les données des axes restent inchangées.

Données pouvant être conservées en dehors du tableau

J'ai mentionné que vous pouvez contenir Coordinate en plus des données de tableau. Les données de coordonnées peuvent être classées comme suit.

Dimension coordinate

Cela a été introduit en tant que Coordinate dans la section précédente. Il porte le même nom que le nom de la dimension («time», «lon», «lat» dans les données de la section précédente).

Non-dimension coordinate

Vous pouvez également contenir des coordonnées avec des noms différents de la dimension. Je vais vous expliquer dans quelles situations il est nécessaire et pratique, mais au fond, il faut comprendre qu'il est ** conservé mais non calculé **.

Scalar coordinate

Vous pouvez contenir un scalaire qui peut être un axe de coordonnées. Il s'agit également d'un type de données fréquemment utilisé, mais pour le moment, vous pouvez le considérer comme une coordonnée sans dimension, qui n'est pas calculée mais qui est conservée.

Méthode de base

Il a les méthodes de base trouvées dans np.ndarray.

autres,

Chose sur le type

Conversion de DataArray en np.ndarray.

Le DataArray contient np.ndarray. (En complément, vous pouvez en fait stocker des tableaux multidimensionnels autres que np.ndarray.)

Vous pouvez accéder au tableau de contenu en utilisant .data. (Notez que .values convertit tout objet dont le contenu est en np.ndarray et le renvoie.)

type(data.data), data.shape
(numpy.ndarray, (2920, 25, 53))

Extraire les informations de l'axe des coordonnées

Les objets d'axe de coordonnées peuvent être récupérés en passant le nom de l'axe dans «[]» comme un dictionnaire.

data['lat']
<xarray.DataArray 'lat' (lat: 25)>
array([75. , 72.5, 70. , 67.5, 65. , 62.5, 60. , 57.5, 55. , 52.5, 50. , 47.5,
       45. , 42.5, 40. , 37.5, 35. , 32.5, 30. , 27.5, 25. , 22.5, 20. , 17.5,
       15. ], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
Attributes:
    standard_name:  latitude
    long_name:      Latitude
    units:          degrees_north
    axis:           Y

Celui récupéré est également un objet DataArray.

xarray.Dataset

Un Dataset est un objet qui est une collection de plusieurs DataArrays.

data = xr.tutorial.open_dataset('air_temperature')
data['mean_temperature'] = data['air'].mean('time')
data
<xarray.Dataset>
Dimensions:           (lat: 25, lon: 53, time: 2920)
Coordinates:
  * lat               (lat) float32 75.0 72.5 70.0 67.5 ... 22.5 20.0 17.5 15.0
  * lon               (lon) float32 200.0 202.5 205.0 ... 325.0 327.5 330.0
  * time              (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Data variables:
    air               (time, lat, lon) float32 241.2 242.5 ... 296.19 295.69
    mean_temperature  (lat, lon) float32 260.37564 260.1826 ... 297.30502
Attributes:
    Conventions:  COARDS
    title:        4x daily NMC reanalysis (1948)
    description:  Data is from NMC initialized reanalysis\n(4x/day).  These a...
    platform:     Model
    references:   http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly...

Dans cet exemple, data contient deux DataArrays, air et mean_temperature. Ils sont répertoriés dans la section intitulée «Variables de données».

est.

Les tableaux de données que vous détenez peuvent partager des axes. Par conséquent, il est possible d'indexer plusieurs DataArrays en même temps, comme indiqué ci-dessous.

data.sel(lat=70, method='nearest')
<xarray.Dataset>
Dimensions:           (lon: 53, time: 2920)
Coordinates:
    lat               float32 70.0
  * lon               (lon) float32 200.0 202.5 205.0 ... 325.0 327.5 330.0
  * time              (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Data variables:
    air               (time, lon) float32 250.0 249.79999 ... 242.59 246.29
    mean_temperature  (lon) float32 264.7681 264.3271 ... 253.58247 257.71475
Attributes:
    Conventions:  COARDS
    title:        4x daily NMC reanalysis (1948)
    description:  Data is from NMC initialized reanalysis\n(4x/day).  These a...
    platform:     Model
    references:   http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly...

Le DataArray conservé peut avoir différentes dimensions. Le fait que les axes soient communs ou non est déterminé par le nom de l'axe. Par conséquent, il n'est pas possible d'avoir des axes de coordonnées avec le même nom, bien qu'ils ne soient pas courants.

Pour obtenir un DataArray à partir d'un Dataset, transmettez simplement le nom du DataArray en tant que clé, comme dans un dictionnaire.

da = data['air']
da
<xarray.DataArray 'air' (time: 2920, lat: 25, lon: 53)>
array([[[241.2    , 242.5    , ..., 235.5    , 238.59999],
        [243.79999, 244.5    , ..., 235.29999, 239.29999],
        ...,
        [295.9    , 296.19998, ..., 295.9    , 295.19998],
        [296.29   , 296.79   , ..., 296.79   , 296.6    ]],

       [[242.09999, 242.7    , ..., 233.59999, 235.79999],
        [243.59999, 244.09999, ..., 232.5    , 235.7    ],
        ...,
        [296.19998, 296.69998, ..., 295.5    , 295.1    ],
        [296.29   , 297.19998, ..., 296.4    , 296.6    ]],

       ...,

       [[245.79   , 244.79   , ..., 243.98999, 244.79   ],
        [249.89   , 249.29   , ..., 242.48999, 244.29   ],
        ...,
        [296.29   , 297.19   , ..., 295.09   , 294.38998],
        [297.79   , 298.38998, ..., 295.49   , 295.19   ]],

       [[245.09   , 244.29   , ..., 241.48999, 241.79   ],
        [249.89   , 249.29   , ..., 240.29   , 241.68999],
        ...,
        [296.09   , 296.88998, ..., 295.69   , 295.19   ],
        [297.69   , 298.09   , ..., 296.19   , 295.69   ]]], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Instanciation DataArray

Cette section décrit comment créer un objet DataArray.

Le moyen le plus simple est de donner np.ndarray et le nom de chaque axe.

xr.DataArray(np.arange(12).reshape(3, 4), dims=['x', 'y'])
<xarray.DataArray (x: 3, y: 4)>
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])
Dimensions without coordinates: x, y

L'exemple ci-dessus crée un DataArray avec des éléments 3x4. Le nom de la première dimension est «x» et la deuxième dimension est «y».

Pour donner les axes, spécifiez le mot-clé coords comme dictionnaire.

xr.DataArray(np.arange(12).reshape(3, 4), dims=['x', 'y'], 
             coords={'x': [0, 1, 2], 'y': [0.1, 0.2, 0.3, 0.4]})
<xarray.DataArray (x: 3, y: 4)>
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])
Coordinates:
  * x        (x) int64 0 1 2
  * y        (y) float64 0.1 0.2 0.3 0.4

Pour inclure une coordonnée sans dimension, vous devez spécifier le nom de l'axe dont elle dépend. Pour ce faire, appuyez sur l'argument de type dictionnaire (à droite des deux points). Le premier élément est le nom de l'axe dont il dépend et le second est le corps du tableau.

xr.DataArray(np.arange(12).reshape(3, 4), dims=['x', 'y'], 
             coords={'x': [0, 1, 2], 'y': [0.1, 0.2, 0.3, 0.4], 
                     'z': ('x', ['a', 'b', 'c'])})
<xarray.DataArray (x: 3, y: 4)>
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])
Coordinates:
  * x        (x) int64 0 1 2
  * y        (y) float64 0.1 0.2 0.3 0.4
    z        (x) <U1 'a' 'b' 'c'

Où la coordonnée non dimensionnelle «z» est définie sur l'axe «x».

Passez également la coordonnée scalaire à coords.

xr.DataArray(np.arange(12).reshape(3, 4), dims=['x', 'y'], 
             coords={'x': [0, 1, 2], 'y': [0.1, 0.2, 0.3, 0.4], 'scalar': 3})
<xarray.DataArray (x: 3, y: 4)>
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])
Coordinates:
  * x        (x) int64 0 1 2
  * y        (y) float64 0.1 0.2 0.3 0.4
    scalar   int64 3

Vous pouvez voir que la coordonnée scalaire n'appartient à aucune dimension.

En plus des informations sur l'axe, vous pouvez également transmettre le nom de l'attribut ou DataArray. Cependant, les attributs sont perdus au milieu de l'opération, il peut donc être préférable de placer les éléments que vous souhaitez enregistrer dans un fichier à la fois.

J'expliquerai cela dans le chapitre sur la structure des données pour le stockage, qui sera expliqué plus tard.

Indexation .isel, .sel, .interp, .reindex

L'indexation est l'opération la plus élémentaire et essentielle. Il n'est pas exagéré de dire que le développement de xarray a commencé à le simplifier.

xarray permet l'indexation basée sur la position et l'indexation basée sur les axes.

Indexation basée sur la position .isel

Ceci est similaire à l'indexation pour les tableaux courants tels que np.ndarray. Donnez un entier comme argument entre crochets, comme data [i, j, k]. L'argument indique la position de l'élément dans le ** tableau **.

da[0, :4, 3]
<xarray.DataArray 'air' (lat: 4)>
array([244.     , 244.2    , 247.5    , 266.69998], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5
    lon      float32 207.5
    time     datetime64[ns] 2013-01-01
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Avec la notation ci-dessus, vous devez vous souvenir de la position de chaque dimension (par exemple, l'axe «temps» est le tout premier axe).

Il est étonnamment difficile de s'en souvenir pour chaque analyse. Si vous l'écrivez ainsi, vous n'avez qu'à vous souvenir du nom de l'axe.

da.isel(time=0, lat=slice(0, 4), lon=3)
<xarray.DataArray 'air' (lat: 4)>
array([244.     , 244.2    , 247.5    , 266.69998], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5
    lon      float32 207.5
    time     datetime64[ns] 2013-01-01
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

De cette façon, mettez le nom de l'axe au format mot-clé dans l'argument de la méthode .isel (), et donnez l'index correspondant après l'égal.

Notez que vous devez explicitement utiliser des tranches (slice (0, 4)) lors de la spécification d'une plage.

Indexation basée sur les coordonnées .sel

Dans l'analyse réelle des données, la plage d'intérêt est souvent spécifiée par des coordonnées. Ceci peut être réalisé avec xarray en utilisant la méthode .sel.

da.sel(time='2013-01-01', lat=slice(75, 60))
<xarray.DataArray 'air' (time: 4, lat: 7, lon: 53)>
array([[[241.2    , 242.5    , ..., 235.5    , 238.59999],
        [243.79999, 244.5    , ..., 235.29999, 239.29999],
        ...,
        [272.1    , 270.9    , ..., 275.4    , 274.19998],
        [273.69998, 273.6    , ..., 274.19998, 275.1    ]],

       [[242.09999, 242.7    , ..., 233.59999, 235.79999],
        [243.59999, 244.09999, ..., 232.5    , 235.7    ],
        ...,
        [269.19998, 268.5    , ..., 275.5    , 274.69998],
        [272.1    , 272.69998, ..., 275.79   , 276.19998]],

       [[242.29999, 242.2    , ..., 236.09999, 238.7    ],
        [244.59999, 244.39   , ..., 232.     , 235.7    ],
        ...,
        [273.     , 273.5    , ..., 275.29   , 274.29   ],
        [275.5    , 275.9    , ..., 277.4    , 277.6    ]],

       [[241.89   , 241.79999, ..., 235.5    , 237.59999],
        [246.29999, 245.29999, ..., 231.5    , 234.5    ],
        ...,
        [273.29   , 272.6    , ..., 277.6    , 276.9    ],
        [274.1    , 274.     , ..., 279.1    , 279.9    ]]], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 62.5 60.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2013-01-01T18:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Comme pour .isel, donnez au mot-clé le nom de l'axe et la valeur de coordonnée souhaitée après l'égal. Vous pouvez également utiliser l'objet slice. Dans ce cas, vous pouvez mettre un nombre réel dans l'argument de slice.

Cependant, il y a quelques réserves à la méthode «.sel».

Lorsque vous souhaitez utiliser des valeurs de coordonnées approximatives

La méthode «.sel» a les options «méthode», «tolérance».

Je pense que la "méthode" la plus couramment utilisée est "la plus proche". Il indexera en utilisant la valeur de coordonnée la plus proche.

da.sel(lat=76, method='nearest')
<xarray.DataArray 'air' (time: 2920, lon: 53)>
array([[241.2    , 242.5    , 243.5    , ..., 232.79999, 235.5    , 238.59999],
       [242.09999, 242.7    , 243.09999, ..., 232.     , 233.59999, 235.79999],
       [242.29999, 242.2    , 242.29999, ..., 234.29999, 236.09999, 238.7    ],
       ...,
       [243.48999, 242.98999, 242.09   , ..., 244.18999, 244.48999, 244.89   ],
       [245.79   , 244.79   , 243.48999, ..., 243.29   , 243.98999, 244.79   ],
       [245.09   , 244.29   , 243.29   , ..., 241.68999, 241.48999, 241.79   ]],
      dtype=float32)
Coordinates:
    lat      float32 75.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Lorsque l'indexation basée sur les axes ne fonctionne pas

Les axes de coordonnées ont un certain degré de liberté, vous ne pourrez donc peut-être pas les indexer correctement. Par exemple, les cas suivants sont applicables.

Lorsque les axes de coordonnées n'augmentent pas ou ne diminuent pas de manière monotone

Vous pouvez trier par cet axe avec la méthode sortby, comme da.sortby ('lat').

Lorsque les axes contiennent Nan

En utilisant la méthode isnull qui renvoie des coordonnées qui ne sont pas Nan, commeda.isel (lat = ~ da ['lat']. Isnull ()), vous pouvez créer un tableau qui ignore ces valeurs de coordonnées. Masu

Lorsqu'il y a des valeurs qui se chevauchent dans les axes

Avec da.isel (lat = np.unique (da ['lat'], return_index = True) [1]), vous pouvez créer un tableau qui saute la partie correspondant aux coordonnées dupliquées.

Interpolation .interp

L'interpolation est possible avec la même syntaxe que l'indexation.

da.interp(lat=74)
<xarray.DataArray 'air' (time: 2920, lon: 53)>
array([[242.23999321, 243.29999998, 243.9799988 , ..., 232.79998779,
        235.41999512, 238.87998962],
       [242.69999081, 243.25999451, 243.53999329, ..., 231.60000001,
        233.1599945 , 235.75999145],
       [243.21998903, 243.07599789, 242.97999266, ..., 232.69998783,
        234.45999444, 237.49999702],
       ...,
       [245.72999275, 245.38999009, 244.68999648, ..., 242.74999088,
        243.20999146, 244.00999451],
       [247.42999566, 246.58999336, 245.48999023, ..., 242.4899933 ,
        243.38999027, 244.58999329],
       [247.00999749, 246.28999329, 245.32999587, ..., 240.84999084,
        241.00999144, 241.74999085]])
Coordinates:
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
    lat      int64 74
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Indexation conditionnelle xr.where

Il y a des moments où vous souhaitez sélectionner uniquement des données répondant à des conditions complexes. Par exemple, dans l'exemple ci-dessus

Ou quelque chose.

Lorsque le résultat de la sélection est unidimensionnel, par exemple, lors de la sélection d'un certain jour, il peut être réalisé en passant un tableau unidimensionnel de type Bool à la méthode ".isel".

import datetime

#Choisissez uniquement les week-ends
is_weekend = da['time'].dt.dayofweek.isin([5, 6])
is_weekend
<xarray.DataArray 'dayofweek' (time: 2920)>
array([False, False, False, ..., False, False, False])
Coordinates:
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
#type booléen.Passez-le à la méthode isel.
da.isel(time=is_weekend)
<xarray.DataArray 'air' (time: 832, lat: 25, lon: 53)>
array([[[248.59999, 246.89   , ..., 245.7    , 246.7    ],
        [256.69998, 254.59999, ..., 248.     , 250.09999],
        ...,
        [296.69998, 296.     , ..., 295.     , 294.9    ],
        [297.     , 296.9    , ..., 296.79   , 296.69998]],

       [[245.09999, 243.59999, ..., 249.89   , 251.39   ],
        [254.     , 251.79999, ..., 247.79999, 250.89   ],
        ...,
        [296.79   , 296.19998, ..., 294.6    , 294.69998],
        [297.19998, 296.69998, ..., 296.29   , 296.6    ]],

       ...,

       [[242.39   , 241.79999, ..., 247.59999, 247.29999],
        [247.18999, 246.09999, ..., 253.29999, 254.29999],
        ...,
        [294.79   , 295.29   , ..., 297.69998, 297.6    ],
        [296.38998, 297.19998, ..., 298.19998, 298.29   ]],

       [[249.18999, 248.     , ..., 244.09999, 242.5    ],
        [254.5    , 253.     , ..., 251.     , 250.59999],
        ...,
        [294.88998, 295.19998, ..., 298.     , 297.29   ],
        [296.1    , 296.88998, ..., 298.5    , 298.29   ]]], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-05 ... 2014-12-28T18:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Lorsque vous souhaitez effectuer une sélection basée sur les données d'un tableau multidimensionnel, vous ne pouvez pas définir la forme du tableau après la sélection, c'est donc une bonne idée de remplacer les éléments non sélectionnés par Nan. Dans ce cas, vous pouvez utiliser la fonction where.

da.where(da > 270)
<xarray.DataArray 'air' (time: 2920, lat: 25, lon: 53)>
array([[[      nan,       nan,       nan, ...,       nan,       nan,
               nan],
        [      nan,       nan,       nan, ...,       nan,       nan,
               nan],
        [      nan,       nan,       nan, ...,       nan,       nan,
               nan],
        ...,
        [296.6    , 296.19998, 296.4    , ..., 295.4    , 295.1    ,
         294.69998],
        [295.9    , 296.19998, 296.79   , ..., 295.9    , 295.9    ,
         295.19998],
        [296.29   , 296.79   , 297.1    , ..., 296.9    , 296.79   ,
         296.6    ]],

       [[      nan,       nan,       nan, ...,       nan,       nan,
               nan],
        [      nan,       nan,       nan, ...,       nan,       nan,
               nan],
        [      nan,       nan,       nan, ...,       nan,       nan,
               nan],
...
        [293.69   , 293.88998, 295.38998, ..., 295.09   , 294.69   ,
         294.29   ],
        [296.29   , 297.19   , 297.59   , ..., 295.29   , 295.09   ,
         294.38998],
        [297.79   , 298.38998, 298.49   , ..., 295.69   , 295.49   ,
         295.19   ]],

       [[      nan,       nan,       nan, ...,       nan,       nan,
               nan],
        [      nan,       nan,       nan, ...,       nan,       nan,
               nan],
        [      nan,       nan,       nan, ...,       nan,       nan,
               nan],
        ...,
        [293.79   , 293.69   , 295.09   , ..., 295.29   , 295.09   ,
         294.69   ],
        [296.09   , 296.88998, 297.19   , ..., 295.69   , 295.69   ,
         295.19   ],
        [297.69   , 298.09   , 298.09   , ..., 296.49   , 296.19   ,
         295.69   ]]], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Ici, les éléments qui ne s'appliquent pas à «da> 270» sont remplacés par «np.nan».

Indexation évolutive

Jusqu'à présent, nous avons introduit l'indexation de base, mais xarray permet également une indexation avancée.

Comme le montre cette figure, vous pouvez découper un autre tableau du tableau multidimensionnel. Dans le document officiel, cela s'appelle l'indexation avancée.

image.png

Pour ce faire, définissez le tableau passé comme argument en tant que DataArray. À ce moment, définissez les dimensions de ces tableaux sur les dimensions du tableau nouvellement créé.

Selon la figure ci-dessus, la nouvelle dimension est «z», définissez donc un DataArray sur «z».

lon_new = xr.DataArray([60, 61, 62], dims=['z'], coords={'z': ['a', 'b', 'c']})
lat_new = xr.DataArray([16, 46, 76], dims=['z'], coords={'z': ['a', 'b', 'c']})

Si vous les passez comme arguments à la méthode .sel,(lon, lat)sera [(60, 16), (61, 46), (62, 76)] (le plus proche de) Choisira.

da.sel(lon=lon_new, lat=lat_new, method='nearest')
<xarray.DataArray 'air' (time: 2920, z: 3)>
array([[296.29   , 280.     , 241.2    ],
       [296.29   , 279.19998, 242.09999],
       [296.4    , 278.6    , 242.29999],
       ...,
       [298.19   , 279.99   , 243.48999],
       [297.79   , 279.69   , 245.79   ],
       [297.69   , 279.79   , 245.09   ]], dtype=float32)
Coordinates:
    lat      (z) float32 15.0 45.0 75.0
    lon      (z) float32 200.0 200.0 200.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
  * z        (z) <U1 'a' 'b' 'c'
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Vous pouvez voir que la séquence nouvellement obtenue ne dépend plus des axes lon, lat, mais du z.

Une méthode similaire peut également être utilisée pour l'interpolation.

da.interp(lon=lon_new, lat=lat_new)
<xarray.DataArray 'air' (time: 2920, z: 3)>
array([[nan, nan, nan],
       [nan, nan, nan],
       [nan, nan, nan],
       ...,
       [nan, nan, nan],
       [nan, nan, nan],
       [nan, nan, nan]])
Coordinates:
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
    lon      (z) int64 60 61 62
    lat      (z) int64 16 46 76
  * z        (z) <U1 'a' 'b' 'c'
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Combiner les données xr.concat

Vous souhaitez souvent combiner plusieurs objets xarray. Par exemple, lorsque vous souhaitez analyser des données expérimentales et des données de simulation obtenues avec divers paramètres ensemble.

C'est l'équivalent de «np.concatenate», «np.stack». Je pense que c'est facile à comprendre si vous y voyez l'opération inverse de l'indexation.

La syntaxe est

xr.concat([data1, data2, ..., ], dim='dim')

est. Passez plusieurs DataArrays (ou Datasets) comme premier argument. Le deuxième argument est la direction de la jointure. En fonction de la direction de jonction, il peut être divisé en opérations de concaténation et d'empilage.

Spécifier un axe existant (comportement de concaténation)

L'opération de base de la fonction xr.concat est de connecter plusieurs DataArrays le long de l'axe existant.

Par exemple, les deux DataArray suivants

da0 = da.isel(time=slice(0, 10))
da1 = da.isel(time=slice(10, None))
da0
<xarray.DataArray 'air' (time: 10, lat: 25, lon: 53)>
array([[[241.2    , 242.5    , ..., 235.5    , 238.59999],
        [243.79999, 244.5    , ..., 235.29999, 239.29999],
        ...,
        [295.9    , 296.19998, ..., 295.9    , 295.19998],
        [296.29   , 296.79   , ..., 296.79   , 296.6    ]],

       [[242.09999, 242.7    , ..., 233.59999, 235.79999],
        [243.59999, 244.09999, ..., 232.5    , 235.7    ],
        ...,
        [296.19998, 296.69998, ..., 295.5    , 295.1    ],
        [296.29   , 297.19998, ..., 296.4    , 296.6    ]],

       ...,

       [[244.79999, 244.39   , ..., 242.7    , 244.79999],
        [246.7    , 247.09999, ..., 237.79999, 240.2    ],
        ...,
        [297.79   , 297.19998, ..., 296.4    , 295.29   ],
        [297.9    , 297.69998, ..., 297.19998, 297.     ]],

       [[243.89   , 243.79999, ..., 240.29999, 242.59999],
        [245.5    , 245.79999, ..., 236.59999, 239.     ],
        ...,
        [297.6    , 297.     , ..., 295.6    , 295.     ],
        [298.1    , 297.69998, ..., 296.79   , 297.1    ]]], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2013-01-03T06:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Vous pouvez restaurer les données d'origine en vous connectant dans le sens «temps».

xr.concat([da0, da1], dim='time')
<xarray.DataArray 'air' (time: 2920, lat: 25, lon: 53)>
array([[[241.2    , 242.5    , 243.5    , ..., 232.79999, 235.5    ,
         238.59999],
        [243.79999, 244.5    , 244.7    , ..., 232.79999, 235.29999,
         239.29999],
        [250.     , 249.79999, 248.89   , ..., 233.2    , 236.39   ,
         241.7    ],
        ...,
        [296.6    , 296.19998, 296.4    , ..., 295.4    , 295.1    ,
         294.69998],
        [295.9    , 296.19998, 296.79   , ..., 295.9    , 295.9    ,
         295.19998],
        [296.29   , 296.79   , 297.1    , ..., 296.9    , 296.79   ,
         296.6    ]],

       [[242.09999, 242.7    , 243.09999, ..., 232.     , 233.59999,
         235.79999],
        [243.59999, 244.09999, 244.2    , ..., 231.     , 232.5    ,
         235.7    ],
        [253.2    , 252.89   , 252.09999, ..., 230.79999, 233.39   ,
         238.5    ],
...
        [293.69   , 293.88998, 295.38998, ..., 295.09   , 294.69   ,
         294.29   ],
        [296.29   , 297.19   , 297.59   , ..., 295.29   , 295.09   ,
         294.38998],
        [297.79   , 298.38998, 298.49   , ..., 295.69   , 295.49   ,
         295.19   ]],

       [[245.09   , 244.29   , 243.29   , ..., 241.68999, 241.48999,
         241.79   ],
        [249.89   , 249.29   , 248.39   , ..., 239.59   , 240.29   ,
         241.68999],
        [262.99   , 262.19   , 261.38998, ..., 239.89   , 242.59   ,
         246.29   ],
        ...,
        [293.79   , 293.69   , 295.09   , ..., 295.29   , 295.09   ,
         294.69   ],
        [296.09   , 296.88998, 297.19   , ..., 295.69   , 295.69   ,
         295.19   ],
        [297.69   , 298.09   , 298.09   , ..., 296.49   , 296.19   ,
         295.69   ]]], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Notez également que l'axe «temps» a également été restauré.

Spécifier un nouvel axe (opération de pile)

Par exemple, supposons que vous ayez deux DataArrays comme celui-ci:

da0 = da.isel(time=0)
da1 = da.isel(time=1)
da0
<xarray.DataArray 'air' (lat: 25, lon: 53)>
array([[241.2    , 242.5    , 243.5    , ..., 232.79999, 235.5    , 238.59999],
       [243.79999, 244.5    , 244.7    , ..., 232.79999, 235.29999, 239.29999],
       [250.     , 249.79999, 248.89   , ..., 233.2    , 236.39   , 241.7    ],
       ...,
       [296.6    , 296.19998, 296.4    , ..., 295.4    , 295.1    , 294.69998],
       [295.9    , 296.19998, 296.79   , ..., 295.9    , 295.9    , 295.19998],
       [296.29   , 296.79   , 297.1    , ..., 296.9    , 296.79   , 296.6    ]],
      dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
    time     datetime64[ns] 2013-01-01
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Pour combiner cela dans la direction de l'axe «temps»:

xr.concat([da0, da1], dim='time')
<xarray.DataArray 'air' (time: 2, lat: 25, lon: 53)>
array([[[241.2    , 242.5    , 243.5    , ..., 232.79999, 235.5    ,
         238.59999],
        [243.79999, 244.5    , 244.7    , ..., 232.79999, 235.29999,
         239.29999],
        [250.     , 249.79999, 248.89   , ..., 233.2    , 236.39   ,
         241.7    ],
        ...,
        [296.6    , 296.19998, 296.4    , ..., 295.4    , 295.1    ,
         294.69998],
        [295.9    , 296.19998, 296.79   , ..., 295.9    , 295.9    ,
         295.19998],
        [296.29   , 296.79   , 297.1    , ..., 296.9    , 296.79   ,
         296.6    ]],

       [[242.09999, 242.7    , 243.09999, ..., 232.     , 233.59999,
         235.79999],
        [243.59999, 244.09999, 244.2    , ..., 231.     , 232.5    ,
         235.7    ],
        [253.2    , 252.89   , 252.09999, ..., 230.79999, 233.39   ,
         238.5    ],
        ...,
        [296.4    , 295.9    , 296.19998, ..., 295.4    , 295.1    ,
         294.79   ],
        [296.19998, 296.69998, 296.79   , ..., 295.6    , 295.5    ,
         295.1    ],
        [296.29   , 297.19998, 297.4    , ..., 296.4    , 296.4    ,
         296.6    ]]], dtype=float32)
Coordinates:
  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 2013-01-01T06:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

De cette façon, s'il y a une coordonnée scalaire appelée «temps» dans «da0», «da1», connectez-la dans cette direction. Vous pouvez faire de l'axe temps une nouvelle coordonnée de dimension (en spécifiant un nom dans le mot-clé dim, tel que dim = 'time').

Vous pouvez également créer un nouvel axe de coordonnées pour ce nom en donnant un DataArray au mot-clé dim.

Fichier IO

Un autre avantage de xarray est qu'il peut stocker des données telles que DataArray et Dataset de manière cohérente et autonome. Il prend en charge plusieurs formats de fichiers, mais le plus simple à utiliser est probablement netCDF.

Wikipedia a également une explication.

NetCDF est un format binaire indépendant du modèle de l'ordinateur (indépendant du modèle), peut lire et écrire des données sous forme de tableau (orienté tableau) et peut stocker des explications sur les données en plus des données (auto-descriptives). ). NetCDF est une norme internationale de l'Open Geospatial Consortium, un consortium international qui développe des normes d'information géographique ouvertes.

Actuellement, les versions 3 et 4 sont en maintenance, mais il est plus sûr d'utiliser la version 4. La version 4 est conforme à la norme HDF5, de sorte que même un lecteur HDF5 ordinaire peut lire le contenu.

L'extension officielle est «.nc».

Enregistrer et charger des fichiers au format netCDF

Le package netcdf4 est requis pour enregistrer au format netCDF.

pip install netcdf4

Ou

conda install netcdf4

Lançons et installons-le sur votre système.

Pour enregistrer l'objet xarray au format netCDF, exécutez la méthode .to_netcdf.

data.to_netcdf('test.nc')
/home/keisukefujii/miniconda3/envs/xarray/lib/python3.7/site-packages/ipykernel_launcher.py:1: SerializationWarning: saving variable air with floating point data as an integer dtype without any _FillValue to use for NaNs
  """Entry point for launching an IPython kernel.

Vous devriez maintenant avoir test.nc enregistré dans votre chemin actuel.

Inversement, pour charger le fichier enregistré, faites xr.load_dataset.

xr.load_dataset('test.nc')
<xarray.Dataset>
Dimensions:           (lat: 25, lon: 53, time: 2920)
Coordinates:
  * lat               (lat) float32 75.0 72.5 70.0 67.5 ... 22.5 20.0 17.5 15.0
  * lon               (lon) float32 200.0 202.5 205.0 ... 325.0 327.5 330.0
  * time              (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Data variables:
    air               (time, lat, lon) float32 241.2 242.5 ... 296.19 295.69
    mean_temperature  (lat, lon) float32 260.37564 260.1826 ... 297.30502
Attributes:
    Conventions:  COARDS
    title:        4x daily NMC reanalysis (1948)
    description:  Data is from NMC initialized reanalysis\n(4x/day).  These a...
    platform:     Model
    references:   http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly...

Comme vous pouvez le voir, non seulement la partie du tableau, mais aussi le nom de l'axe, les valeurs de coordonnées et d'autres informations dans les attributs sont stockés.

Une note sur scipy.io.netcdf

Le xarray .to_netcdf utilise la bibliothèque netcdf4 par défaut, mais utilise scipy.io.netcdf s'il n'est pas installé sur votre système. Cependant, scipy.io.netcdf prend en charge la version 3, qui est incompatible avec la version 4 et peut prêter à confusion.

N'oubliez pas d'installer netcdf4 sur votre système. Vous pouvez également spécifier le package à utiliser explicitement, tel que .to_netcdf ('test.nc', engine = 'netcdf4').

Baie de mémoire insuffisante

Les fichiers au format NetCDF peuvent être lus de manière aléatoire. Par conséquent, si vous souhaitez utiliser uniquement les données d'une certaine partie du fichier, il est efficace de lire toutes les données lorsque cela est nécessaire au lieu de les lire d'abord à partir du disque dur. Surtout lorsqu'il s'agit de fichiers trop volumineux pour être stockés en mémoire, un tel traitement devient essentiel. C'est ce qu'on appelle le traitement hors mémoire.

Le traitement du manque de mémoire dans xarray est positionné comme un hier assez important, et diverses implémentations sont en cours. Je présenterai les détails dans un autre article, et je ne mentionnerai ici que les contenus les plus élémentaires.

Le traitement de manque de mémoire est effectué en faisant open_dataset au lieu de load_dataset lors du chargement.

unloaded_data = xr.open_dataset('test.nc')
unloaded_data
<xarray.Dataset>
Dimensions:           (lat: 25, lon: 53, time: 2920)
Coordinates:
  * lat               (lat) float32 75.0 72.5 70.0 67.5 ... 22.5 20.0 17.5 15.0
  * lon               (lon) float32 200.0 202.5 205.0 ... 325.0 327.5 330.0
  * time              (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Data variables:
    air               (time, lat, lon) float32 ...
    mean_temperature  (lat, lon) float32 ...
Attributes:
    Conventions:  COARDS
    title:        4x daily NMC reanalysis (1948)
    description:  Data is from NMC initialized reanalysis\n(4x/day).  These a...
    platform:     Model
    references:   http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly...

Si vous regardez "Variables de données", par exemple, vous pouvez voir qu'il n'y a pas de nombre dans "air", c'est juste "...".

Cela indique qu'aucune donnée n'a été lue. Cependant, comme les axes de coordonnées sont toujours lus, le travail d'indexation ci-dessus peut être effectué tel quel.

unloaded_data['air'].sel(lat=60, method='nearest')
<xarray.DataArray 'air' (time: 2920, lon: 53)>
[154760 values with dtype=float32]
Coordinates:
    lat      float32 60.0
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Attributes:
    long_name:     4xDaily Air temperature at sigma level 995
    units:         degK
    precision:     2
    GRIB_id:       11
    GRIB_name:     TMP
    var_desc:      Air temperature
    dataset:       NMC Reanalysis
    level_desc:    Surface
    statistic:     Individual Obs
    parent_stat:   Other
    actual_range:  [185.16 322.1 ]

Lorsque je sélectionne une variable ou que je l'indexe, elle ne se charge toujours pas. Il sera chargé en mémoire pour la première fois lorsqu'il sera calculé ou lorsque .data ou .values sera appelé.

Vous pouvez également utiliser la méthode .compute () pour le charger explicitement.

unloaded_data.compute()
<xarray.Dataset>
Dimensions:           (lat: 25, lon: 53, time: 2920)
Coordinates:
  * lat               (lat) float32 75.0 72.5 70.0 67.5 ... 22.5 20.0 17.5 15.0
  * lon               (lon) float32 200.0 202.5 205.0 ... 325.0 327.5 330.0
  * time              (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Data variables:
    air               (time, lat, lon) float32 241.2 242.5 ... 296.19 295.69
    mean_temperature  (lat, lon) float32 260.37564 260.1826 ... 297.30502
Attributes:
    Conventions:  COARDS
    title:        4x daily NMC reanalysis (1948)
    description:  Data is from NMC initialized reanalysis\n(4x/day).  These a...
    platform:     Model
    references:   http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly...

Remarque

open_dataset, open_dataarray verrouille le fichier pour une lecture future. Par conséquent, si vous touchez un fichier du même nom dans un autre programme, une erreur se produira.

Par exemple, considérons le cas où le résultat d'un calcul chronophage est généré à chaque fois sous forme de fichier «.nc». Si vous voulez voir la progression et l'ouvrir avec open_dataset, toutes les tentatives ultérieures d'écriture dans ce fichier échoueront. C'est mieux s'il échoue et s'arrête, mais il semble que cela détruit souvent les données existantes. Ce processus de mémoire insuffisante doit être utilisé avec prudence.

Vous pouvez utiliser l'instruction with pour vous assurer que le fichier est fermé et déverrouillé, comme indiqué ci-dessous.

with xr.open_dataset('test.nc') as f:
  print(f)
<xarray.Dataset>
Dimensions:           (lat: 25, lon: 53, time: 2920)
Coordinates:
  * lat               (lat) float32 75.0 72.5 70.0 67.5 ... 22.5 20.0 17.5 15.0
  * lon               (lon) float32 200.0 202.5 205.0 ... 325.0 327.5 330.0
  * time              (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00
Data variables:
    air               (time, lat, lon) float32 ...
    mean_temperature  (lat, lon) float32 ...
Attributes:
    Conventions:  COARDS
    title:        4x daily NMC reanalysis (1948)
    description:  Data is from NMC initialized reanalysis\n(4x/day).  These a...
    platform:     Model
    references:   http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly...

Structure de données recommandée

Comme mentionné ci-dessus, si vous enregistrez les données au format netcdf, vous pouvez analyser les données en douceur avec xarray. Lors de la création de données pour des expériences ou des simulations, il est judicieux de créer une structure de données qui prend en compte l'analyse ultérieure et de l'enregistrer au format netcdf en l'incluant.

Il est recommandé de conserver les données suivantes supprimées.

Dans le cas d'une expérience, je pense qu'elle devrait être sauvegardée avec la structure suivante.

from datetime import datetime
#Supposons qu'il s'agit de données d'image prises avec un appareil photo.
raw_data = np.arange(512 * 2048).reshape(512, 2048)
#Coordonnées des axes x et y de l'image.
x = np.arange(512)
y = np.arange(2048)

#Données de position pour chaque pixel de l'image
r = 300.0
X = r * np.cos(x / 512)[:, np.newaxis] + r * np.sin(y / 2048)
Y = r * np.sin(x / 512)[:, np.newaxis] - r * np.cos(y / 2048)

#Paramètres utilisés dans l'expérience
p0 = 3.0
p1 = 'normal'

#Durée de l'expérience
now = datetime.now()

data = xr.DataArray(
    raw_data, dims=['x', 'y'],
    coords={'x': x, 'y': y, 
            'X': (('x', 'y'), X),
            'Y': (('x', 'y'), Y),
            'p0': p0, 'p1': p1,
            'datetime': now},
    attrs={'camera type': 'flash 4', 
           'about': 'Example of the recommended data structure.'})
data
<xarray.DataArray (x: 512, y: 2048)>
array([[      0,       1,       2, ...,    2045,    2046,    2047],
       [   2048,    2049,    2050, ...,    4093,    4094,    4095],
       [   4096,    4097,    4098, ...,    6141,    6142,    6143],
       ...,
       [1042432, 1042433, 1042434, ..., 1044477, 1044478, 1044479],
       [1044480, 1044481, 1044482, ..., 1046525, 1046526, 1046527],
       [1046528, 1046529, 1046530, ..., 1048573, 1048574, 1048575]])
Coordinates:
  * x         (x) int64 0 1 2 3 4 5 6 7 8 ... 504 505 506 507 508 509 510 511
  * y         (y) int64 0 1 2 3 4 5 6 7 ... 2041 2042 2043 2044 2045 2046 2047
    X         (x, y) float64 300.0 300.1 300.3 300.4 ... 414.7 414.8 414.9 414.9
    Y         (x, y) float64 -300.0 -300.0 -300.0 -300.0 ... 89.66 89.79 89.91
    p0        float64 3.0
    p1        <U6 'normal'
    datetime  datetime64[ns] 2020-11-13T03:31:56.219372
Attributes:
    camera type:  flash 4
    about:        Example of the recommended data structure.

En regroupant les choses qui peuvent changer au cours de l'expérience en coordonnées scalaires, il sera plus facile d'analyser plus tard lorsque vous voudrez connaître les dépendances sur eux.

Autres fonctions / méthodes utiles

Reduction .sum, .mean, etc

Des processus de réduction tels que "np.sum" et "np.mean" sont fournis comme méthodes. Par exemple, np.sum vous permet de spécifier dans quelle direction d'axe ajouter. De même, xarray vous permet de spécifier un axe comme nom d'axe.

data = xr.tutorial.open_dataset('air_temperature')['air']

data.sum('lat')
<xarray.DataArray 'air' (time: 2920, lon: 53)>
array([[6984.9497, 6991.6606, 6991.5303, ..., 6998.77  , 7007.8804,
        7016.5605],
       [6976.4307, 6988.45  , 6993.2407, ..., 6994.3906, 7006.7505,
        7019.941 ],
       [6975.2603, 6982.02  , 6988.77  , ..., 6992.0503, 7004.9404,
        7020.3506],
       ...,
       [6990.7505, 6998.3496, 7013.3496, ..., 6995.05  , 7008.6504,
        7019.4497],
       [6984.95  , 6991.6504, 7007.949 , ..., 6994.15  , 7008.55  ,
        7020.8506],
       [6981.75  , 6983.85  , 6997.0503, ..., 6985.6494, 6999.2495,
        7012.0493]], dtype=float32)
Coordinates:
  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0
  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00

Si vous additionnez dans une certaine direction axiale de cette manière, il n'y aura pas d'axes de coordonnées dans cette direction. Les autres axes restent.

Lors de l'exécution de .sum etc., np.nansum est appelé à la place de np.sum par défaut, et le calcul sautant" Nan "est effectué. Si vous ne voulez pas sauter Nan (si vous voulez utiliser np.sum au lieu de np.nansum), vous devez donner skipna = False.

Dessin .plot

DataArray a une méthode plot, ce qui facilite la visualisation du contenu des données.

data.sel(time='2014-08-30 08:00', method='nearest').plot()
<matplotlib.collections.QuadMesh at 0x7f9a98557310>

output_102_1.png

Il peut être préférable de le considérer comme simple plutôt que de faire un schéma final du document / rapport, mais afin de visualiser et de comprendre les données avec Jupyter, Ipython, etc. et de procéder à l'analyse. Très pratique à utiliser.

Recommended Posts

Analyse des données à l'aide de xarray
Xarray de bibliothèque d'analyse de données multidimensionnelle
Analyse de données à l'aide de pandas python
Recommandation d'analyse des données à l'aide de MessagePack
Bibliothèque d'analyse de données multidimensionnelle xarray Partie 2
Analyse des données Titanic 2
Analyse de données python
Analyse des données Titanic 1
Analyse des données Titanic 3
Création d'une application d'analyse de données à l'aide de Streamlit
Analyse de données avec python 2
Présentation de l'analyse de données python
Filtre médian utilisant xarray (filtre médian)
Nettoyage des données 2 Nettoyage des données à l'aide de DataFrame
Nettoyage des données à l'aide de Python
Modèle d'analyse de données Python
Analyse orthologue à l'aide d'OrthoFinder
Analyse de données avec Python
[Python] [Word] [python-docx] Analyse simple des données de diff en utilisant python
Analyse de Big Data à l'aide du framework de contrôle de flux de données Luigi
[Livre technique] Introduction à l'analyse de données avec Python -1 Chapitre Introduction-
Mon conteneur d'analyse de données python
Python pour l'analyse des données Chapitre 4
Analyse morphologique japonaise avec Janome
Sélectionnez des fonctionnalités avec des données textuelles
[Python] Notes sur l'analyse des données
Notes d'apprentissage sur l'analyse des données Python
Méthode de visualisation de données utilisant matplotlib (1)
Méthode de visualisation de données utilisant matplotlib (2)
Python pour l'analyse des données Chapitre 2
Wrap Analysis part1 (préparation des données)
Conseils et précautions lors de l'analyse des données
Python pour l'analyse des données Chapitre 3
Analyse des données Twitter | Analyse des tendances
[Introduction] Analyse de données satellitaires artificielles à l'aide de Python (environnement Google Colab)
J'ai essayé d'analyser les données scRNA-seq en utilisant l'analyse des données topologiques (TDA)
Première analyse de données satellitaires par Tellus
Obtenir des données Salesforce à l'aide de l'API REST
Acquisition de données à l'aide de l'API googlemap de python
Précautions lors de l'utilisation de l'analyse des traits TextBlob
Méthode de visualisation de données utilisant matplotlib (+ pandas) (5)
Python: analyse des séries chronologiques: prétraitement des données des séries chronologiques
Analyser les données au format CSV à l'aide de SQL
Reconnaissance faciale à l'aide de l'analyse des composants principaux
Obtenez des données Amazon à l'aide de Keep API # 1 Obtenez des données
Méthode de visualisation de données utilisant matplotlib (+ pandas) (3)
Modèle de prétraitement pour l'analyse des données (Python)
Mémo d'acquisition de données à l'aide de l'API Backlog
Expérience de réussite du test d'analyse des données de la version de novembre 2020
Analyse de données pour améliorer POG 3 ~ Analyse de régression ~
Traitement de l'analyse japonaise à l'aide de Janome part1
Binarisation d'images par analyse discriminante linéaire
Analyse des séries chronologiques 3 Prétraitement des données des séries chronologiques
Obtenez des données de Twitter avec Tweepy
Analyse de données à partir de python (visualisation de données 1)
Méthode de visualisation de données utilisant matplotlib (+ pandas) (4)
Analyse de données à partir de python (visualisation de données 2)
Tutoriel de recommandation utilisant l'analyse d'association (concept)
Traitement des données 2 Analyse de divers formats de données