Transpose of a very large matrix with fewer than 1021 rows with Python

Question

I have a matrix of fewer than 1021 rows. I'd like to take its transpose, but it is too big for GNU datamash and awk solutions to fit into memory.

To keep memory overhead low, my thought is that I keep a list of read-only file pointers to the start of each row. I then read out bytes from each pointer until I hit a delimiter (tab or newline).

Because I have fewer than 1021 rows, I won't hit the usual 1024 OS-based file pointer limit.

Once I have read out a field from all file pointers, I write a newline character and start over, reading a field from all file pointers, and again, until no more bytes are available:

#!/usr/bin/env python                                                                                                                                                                                                                               

import sys
import os

try:
  in_fn = sys.argv[1]
except ValueError as ve:
  sys.exit(-1)

def get_size(fn):
  st = os.stat(fn)
  return st.st_size

#                                                                                                                                                                                                                                                   
# 1) Read in file offsets to start of new line                                                                                                                                                                                                      
# 2) Open up a file pointer to that offset                                                                                                                                                                                                          
# 3) Process each file pointer to get a value until a delimiter is hit, then write it as a row of output                                                                                                                                            
#                                                                                                                                                                                                                                                   

size = get_size(in_fn)

fps = []
new_fp = open(in_fn, 'r')
new_fp.seek(0, 0)
fps.append(new_fp)

with open(in_fn, 'r') as f:
  byte = f.read(1)
  while byte:
    if byte == '\n':
      new_offset = f.tell()
      if new_offset < size:
        new_fp = open(in_fn, 'r')
        new_fp.seek(new_offset, 0)
        fps.append(new_fp)
      else:
        break
    byte = f.read(1)

while size > 0:
  for fi, f in enumerate(fps):
    byte = f.read(1)
    size -= 1
    while byte:
      if byte == '\t' or byte == '\n':
        if fi != len(fps) - 1: sys.stdout.write('\t')
        break
      sys.stdout.write('%s' % (byte))
      byte = f.read(1)
      size -= 1
  sys.stdout.write('\n')

for f in fps:
  f.close()

Is there anything I can do to improve the performance of this? Reading and processing a set of file pointers one byte at a time seems quite expensive. However, I need to find newline characters to build offsets and file pointers. Is there a cleverer/faster way (in Python) to find the byte offsets of newlines?

Answer 1

The golden rule of I/O performance is never read a block twice. You're reading each block thousands of times.

I made a thousand-row test file with records of about 2MB each:

tr -dc '[:alnum:]' < /dev/urandom  | fold -w2000000 | head -1000 > test.txt

And modified the fps-building portion of your code thusly:

#!/usr/bin/env python
from prettyprinter import pprint
import sys
import os
import re

try:
  in_fn = sys.argv[1]
except ValueError as ve:
  sys.exit(-1)

sz = os.path.getsize(in_fn)
fps = [ open(in_fn, 'r') ]
bufsz = 100 * 2**20 # MB

with open(in_fn, 'r') as f:
  offset = 0
  buf = f.read(bufsz)
  while len(buf)>0:
    for p in [m.start() for m in re.finditer(r'\n', buf)]:
      new_pos = offset+p+1
      if new_pos < sz:
        new_fp = open(in_fn, 'r')
        new_fp.seek(new_pos, 0)
        fps.append(new_fp)

    offset = f.tell()
    buf = f.read(bufsz)

pprint( list( map( lambda f: f.tell(), fps ) ) )

This runs 65 times faster than the bytewise version on my machine.

The actual transposition will benefit similarly from this treatment. Keep a buffer for each row and append to it with a read when it's empty.

If you want to keep the logic simpler, you can compromise (with a performance penalty) by reading a few kilobytes at a time and rewind-seeking the filehandle back to the first tab character. This would remove the need to keep an array of buffers.