Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Access data from HDF5 - slice/extract data

I am trying to target specific rows within a large matrix contained in an HDF5 file. I have a .txt file that contains the ids of interest that are present in the HDF5 file and wish to output the corresponding data of those rows - all corresponding data are numerals.

I wrote the following code but the output only contains the ids (single column). I need the remaining data attached to those rows (the data in the subsequent columns). Any advice would be much appreciated!

import os 
import h5py


mydir = os.path.expanduser("~/Desktop/alexs-stuff/")

in_file = mydir + "EMP/EMPopen/full_emp_table_hdf5.h5"
wanted_file = mydir + "EMP/greengenes-curto-only.txt"
out_file = mydir + "EMP/emp-curto-only.txt"

wanted = set()

with open(wanted_file) as f:
    for line in f:
        line = line.strip()
        if line != "":
            wanted.add(line)

hdf5_file = h5py.File(in_file, "r")

count = 0

with open(out_file, "w") as h:
    for keys in hdf5_file["observation"]["ids"]:
        if keys in wanted:
            count = count + 1
            h.write(keys + "\n")

print "Converted %i records" % count 

hdf5_file.close()

If it helps, here is the structure of the hdf5 file:

<HDF5 file "full_emp_table_hdf5.h5" (mode r)> (File) /
 sample     /sample (Group) /sample
     metadata     /sample/metadata (Group) /sample/metadata
     group-metadata     /sample/group-metadata (Group) /sample/group-metadata
     ids     /sample/ids (Dataset) /sample/ids     len = (15481,) object
     matrix     /sample/matrix (Group) /sample/matrix
         indices     /sample/matrix/indices (Dataset) /sample/matrix/indices     len = (107439386,) int32
         indptr     /sample/matrix/indptr (Dataset) /sample/matrix/indptr     len = (15482,) int32
         data     /sample/matrix/data (Dataset) /sample/matrix/data     len = (107439386,) float64
 observation     /observation (Group) /observation
     metadata     /observation/metadata (Group) /observation/metadata
         taxonomy     /observation/metadata/taxonomy (Dataset) /observation/metadata/taxonomy     len = (5594412, 7) object
     group-metadata     /observation/group-metadata (Group) /observation/group-metadata
     ids     /observation/ids (Dataset) /observation/ids     len = (5594412,) object
     matrix     /observation/matrix (Group) /observation/matrix
         indices     /observation/matrix/indices (Dataset) /observation/matrix/indices     len = (107439386,) int32
         indptr     /observation/matrix/indptr (Dataset) /observation/matrix/indptr     len = (5594413,) int32
         data     /observation/matrix/data (Dataset) /observation/matrix/data     len = (107439386,) float64

additional information:

type(hdf5_file['observation']['ids'])
>>> <class 'h5py._hl.dataset.Dataset'>


dir(hdf5_file['observation']['ids'])
>>> ['__array__', '__class__', '__delattr__', '__dict__', '__doc__', '__eq__', '__format__', '__getattribute__', '__getitem__', '__hash__', '__init__', '__iter__', '__len__', '__module__', '__ne__', '__new__', '__nonzero__', '__reduce__', '__reduce_ex__', '__repr__', '__setattr__', '__setitem__', '__sizeof__', '__str__', '__subclasshook__', '__weakref__', '_d', '_dcpl', '_e', '_filters', '_id', '_lapl', '_lcpl', '_local', 'astype', 'attrs', 'chunks', 'compression', 'compression_opts', 'dims', 'dtype', 'file', 'fillvalue', 'fletcher32', 'id', 'len', 'maxshape', 'name', 'parent', 'read_direct', 'ref', 'regionref', 'resize', 'scaleoffset', 'shape', 'shuffle', 'size', 'value', 'write_direct']
like image 287
Alex Chase Avatar asked Jul 15 '26 11:07

Alex Chase


1 Answers

I think this hdf5 file was produced by biom-format software, http://biom-format.org, From its pypi package: https://pypi.python.org/packages/source/b/biom-format/biom-format-2.1.2.tar.gz

I can extract a sample hdf5 file: examples/min_sparse_otu_table_hdf5.biom, and process it with:

import numpy as np
from scipy import sparse as sp
import h5py

f = h5py.File('min_sparse_otu_table_hdf5.biom','r')
f.visit(print)  # to see its structure
"""
observation
observation/group-metadata
observation/ids
observation/matrix
observation/matrix/data
observation/matrix/indices
observation/matrix/indptr
observation/metadata
sample
sample/group-metadata
sample/ids
sample/matrix
sample/matrix/data
sample/matrix/indices
sample/matrix/indptr
sample/metadata
"""

Or in more detail:

f.visititems(lambda name,obj:print(name, obj))
"""
...
sample <HDF5 group "/sample" (4 members)>
sample/group-metadata <HDF5 group "/sample/group-metadata" (0 members)>
sample/ids <HDF5 dataset "ids": shape (6,), type "|O4">
sample/matrix <HDF5 group "/sample/matrix" (3 members)>
sample/matrix/data <HDF5 dataset "data": shape (15,), type "<f8">
sample/matrix/indices <HDF5 dataset "indices": shape (15,), type "<i4">
sample/matrix/indptr <HDF5 dataset "indptr": shape (7,), type "<i4">
sample/metadata <HDF5 group "/sample/metadata" (0 members)>
"""

This has a similar structure to your full_emp_table_hdf5.h5

Load arrays with .value (or [:])

sample_ids = f['sample/ids'].value
"""
array([b'Sample1', b'Sample2', b'Sample3', b'Sample4', b'Sample5',
       b'Sample6'], dtype=object)
"""

# create a sparse array from the sample/matrix group
data = f['sample/matrix/data'].value
indices = f['sample/matrix/indices'].value
indptr = f['sample/matrix/indptr'].value
sample = sp.csr_matrix((data,indices,indptr))
# display as a dense array
print(sample.A)

""" 
array([[ 0.,  5.,  0.,  2.,  0.],
       [ 0.,  1.,  0.,  1.,  1.],
       [ 1.,  0.,  1.,  1.,  1.],
       [ 0.,  2.,  4.,  0.,  0.],
       [ 0.,  3.,  0.,  0.,  0.],
       [ 0.,  1.,  2.,  1.,  0.]])
"""

# select a specific row
print(sample[2,:])
"""
  (0, 0)    1.0
  (0, 2)    1.0
  (0, 3)    1.0
  (0, 4)    1.0
"""
print(sample[2,:].A)  # that row in dense format
# array([[ 1.,  0.,  1.,  1.,  1.]])

I can write sample to a text file with:

np.savetxt('min_sample.txt',sample.A, '%.2f')

"""
0.00 5.00 0.00 2.00 0.00
0.00 1.00 0.00 1.00 1.00
1.00 0.00 1.00 1.00 1.00
0.00 2.00 4.00 0.00 0.00
0.00 3.00 0.00 0.00 0.00
0.00 1.00 2.00 1.00 0.00
"""

The biom module probably can do all this and more, but the basics just need the numpy and sparse modules.

like image 180
hpaulj Avatar answered Jul 18 '26 01:07

hpaulj