Проблемы обработки большого набора данных

Решение технических вычислительных задач, требующих обработки и анализа больших объемов данных, ставит высокий спрос на вашу компьютерную систему. Большие наборы данных занимают значительную память во время обработки и могут потребовать многих операций для вычисления решения. Также может потребоваться много времени, чтобы получить доступ к информации из больших файлов данных.

Компьютерные системы, однако, имеют ограниченную память и конечную скорость центральный процессор. Доступные ресурсы варьируются в зависимости от процессора и операционной системы, последняя из которых также потребляет ресурсы. Для примера:

32-разрядные процессоры и операционные системы могут обрабатывать до 2 ^ 32 = 4 294 967 296 = 4 ГБ памяти (также известной как виртуальное адресное пространство). Windows ® XP и Windows ® 2000 выделяют только 2 ГБ этой виртуальной памяти каждому процессу (такому как MATLAB). В UNIX ® виртуальная память, выделенная процессу, конфигурируется системой и обычно составляет около 3 ГБ. Приложение, выполняющее вычисление, такое как MATLAB, может потребовать хранения в дополнение к пользовательской задаче. Основной проблемой при обработке больших объемов данных является то, что требования к памяти программы могут превысить доступные на платформе. Например, MATLAB генерирует ошибку «нехватки памяти», когда требования к данным превышают приблизительно 1,7 ГБ в Windows XP.

Для получения дополнительной информации об управлении памятью и больших наборах данных смотрите Производительность и Память.

На типовой 32-битной машине максимальный размер одного набора данных, с которым можно работать в MATLAB, составляет несколько сотен МБ, или около размера большой хромосомы. Размещение в ОЗУ файлов позволяет MATLAB обойти это ограничение и позволяет вам работать с очень большими наборами данных интуитивно понятным способом.

Целые наборы данных генома

Последние наборы данных генома можно загрузить с веб-сайта Ensembl. Данные предоставляются в нескольких форматах. Они регулярно обновляются по мере появления новой информации о последовательности. Этот пример будет использовать данные ДНК человека, хранящиеся в формате FASTA. Chromosome 1 является (в GRCh37.56 релизе сентября 2009 года) сжатым файлом на 65,6 МБ. После распаковки файла речь идет о 250MB. MATLAB использует 2 байта на символ, поэтому, если вы считаете файл в MATLAB, то это потребует около 500MB памяти.

Этот пример предполагает, что вы уже загрузили и не сжали файл FASTA в свою локальную директорию. Измените имя переменной FASTAfilename при необходимости.

FASTAfilename = 'Homo_sapiens.GRCh37.56.dna.chromosome.1.fa';
fileInfo = dir(which(FASTAfilename))

fileInfo = 

  struct with fields:

       name: 'Homo_sapiens.GRCh37.56.dna.chromosome.1.fa'
     folder: '/mathworks/hub/qe/test_data/Bioinformatics_Toolbox/v000/demoData/biomemorymapdemo'
       date: '01-Feb-2013 11:54:41'
      bytes: 253404851
      isdir: 0
    datenum: 7.3527e+05

Область Отображаемых файлов

Размещение в ОЗУ позволяет MATLAB получить доступ к данным в файле, как будто он находится в памяти. Для доступа к данным можно использовать стандартные операции индексации MATLAB. См. документацию для memmapfile для получения дополнительной информации.

Можно просто сопоставить файл FASTA и получить доступ к данным непосредственно оттуда. Однако файл формата FASTA содержит новые символы линий. The memmapfile функция обрабатывает эти символы так же, как и все другие символы. Удаление их перед размещением в ОЗУ в файле упростит операции индексации. Кроме того, размещение в ОЗУ не работает непосредственно с символом данными, поэтому вам придется рассматривать данные как 8-битные целые числа (класс uint8). Функция nt2int в Bioinformatics Toolbox™ может использоваться, чтобы преобразовать информацию о символах в целочисленные значения. int2nt используется для преобразования назад в символы.

Сначала откройте файл FASTA и извлеките заголовок.

fidIn = fopen(FASTAfilename,'r');
header = fgetl(fidIn)

header =

    '>1 dna:chromosome chromosome:GRCh37:1:1:249250621:1'

Откройте файл для сопоставления памяти.

[fullPath, filename, extension] = fileparts(FASTAfilename);
mmFilename = [filename '.mm']
fidOut = fopen(mmFilename,'w');

mmFilename =

    'Homo_sapiens.GRCh37.56.dna.chromosome.1.mm'

Чтение файла FASTA в блоках 1MB, удаление новых линий, преобразование в uint8 и запись в файл MM.

newLine = sprintf('\n');
blockSize = 2^20;
while ~feof(fidIn)
    % Read in the data
    charData = fread(fidIn,blockSize,'*char')';
    % Remove new lines
    charData = strrep(charData,newLine,'');
    % Convert to integers
    intData = nt2int(charData);
    % Write to the new file
    fwrite(fidOut,intData,'uint8');
end

Закройте файлы.

fclose(fidIn);
fclose(fidOut);

Новый файл имеет примерно тот же размер, что и старый файл, но не содержит новых линий или информации о заголовке.

mmfileInfo = dir(mmFilename)

mmfileInfo = 

  struct with fields:

       name: 'Homo_sapiens.GRCh37.56.dna.chromosome.1.mm'
     folder: '/tmp/BR2021ad_1584584_202060/publish_examples0/tp07af436a/ex57563178'
       date: '26-Jan-2021 19:40:46'
      bytes: 249250621
      isdir: 0
    datenum: 7.3818e+05

Доступ к данным в отображаемом файле памяти

Область команды memmapfile создает объект memmapfile, который сопоставляет новый файл с памятью. В порядок для этого ему нужно знать формат файла. Формат этого файла прост, хотя могут быть сопоставлены намного более сложные форматы.

chr1 = memmapfile(mmFilename, 'format', 'uint8')

chr1 = 

    Filename: '/tmp/BR2021ad_1584584_202060/publish_examples0/tp07af436a/ex57563178/Homo_sapiens.GRCh37.56.dna.chromosome.1.mm'
    Writable: false
      Offset: 0
      Format: 'uint8'
      Repeat: Inf
        Data: 249250621x1 uint8 array

Объект MEMMAPFILE

Объект memmapfile имеет различные свойства. Filename сохраняет полный путь к файлу. Writable указывает, можно ли изменять данные. Обратите внимание, что если вы измените данные, это также изменит исходный файл. Offset позволяет вам задавать пространство, используемое любой информацией о заголовке. Format указывает формат данных. Repeat используется, чтобы задать, сколько блоков (как определяется Format) на карту. Это может быть полезно для ограничения объема памяти, используемой для создания карты памяти. К этим свойствам можно обращаться так же, как к другим данным MATLAB. Для получения дополнительной информации смотрите help memmapfile типа или doc memmapfile.

chr1.Data(1:10)

ans =

  10x1 uint8 column vector

   15
   15
   15
   15
   15
   15
   15
   15
   15
   15

Вы можете получить доступ к любой области данных с помощью операций индексации.

chr1.Data(10000000:10000010)'

ans =

  1x11 uint8 row vector

   1   1   2   2   2   2   3   4   2   4   2

Помните, что нуклеотидная информация была преобразована в целые числа. Можно использовать int2nt для получения информации о последовательности.

int2nt(chr1.Data(10000000:10000010)')

ans =

    'AACCCCGTCTC'

Или используйте seqdisp для отображения последовательности.

seqdisp(chr1.Data(10000000:10001000)')

ans =

  17x71 char array

    '   1  AACCCCGTCT CTACAATAAA TTAAAATATT AGCTGGGCAT GGTGGTGTGT GCTTGTAGTC'
    '  61  CCAGCTACTT GGCGGGCTGA GGTGGGAGAA TCATCCAAGC CTTGGAGGCA GAGGTTGCAG'
    ' 121  TGAGCTGAGA TTGTGACACT GCACTCCAGC CTGGGAGACA GAGTGAGACT CCTACTCAAA'
    ' 181  AAAAAACAAA AAACAAAAAA CAAACCACAA AACTTTCCAG GTAACTTATT AAAACATGTT'
    ' 241  TTTTGTTTGT TTTGAGACAG AGTCTTGCTC TGTCGCCCAG GCTGGAGTGC AGTGGAGCAA'
    ' 301  TCTCAGCTCA CTGCAAGCTC CGCCTCCCGG GTTCACACCA TTCTCCTGCC TCAGCCTCCC'
    ' 361  GAGTAGCTAG GACTATAGGC ACCCGCCACC ACGCCCAGCT TATTTTTTTT GTATTTTTTA'
    ' 421  GTAGAGACGG GGTTTCATCG TGTTAGCCAG GATGGTCTCG ATCTCCTGAC CTCGTGATCC'
    ' 481  GCCCACCTCA GCCTCCCAAA GTGCTGGGAT TACAGGCGTG AGCCACTGCA CCCGGCCTAG'
    ' 541  TTTTTGTATA TTTTTTTTAG TAGAGACAGG GTTTCACCAT GTTAGCCAGG ATGGTCTCAA'
    ' 601  TCTCCTGACC TCGTGATCCG CCCGCCTCGG CCTCCCAAAG TGCTGGGGTT ACAGGCGTGA'
    ' 661  GCCACCGCAC ACAGCATTAA AGCATGTTTT ATTTTCCTAC ACATAATGAA ATCATTACCA'
    ' 721  GATGATTTGA CATGTGTACT TCATTGGAGA GGATTCTTAC AGTATATTCA AAATTAAATA'
    ' 781  TAATGACAAA AAATTACTAC CTAATCTATT AAAATTGGCA TAAGTCATCT ATGATCATTA'
    ' 841  ATGATATGCA AACATAAACA AGTATTATAC CCAGAAGTGT AATTTATTGT AGCTACATCT'
    ' 901  TATGTATAAT AGTTTAGTGG ATTTTTCCTG GAAATTGTCC ATTTTAATTT TTCTCTTAAG'
    ' 961  TCTGTGGAAT TTTCCAGTAA AAGTCAAGGC AAACCCAAGA T                    '

Анализ всей хромосомы

Теперь, когда вы можете легко получить доступ ко всей хромосоме, можно проанализировать данные. Этот пример показывает один из способов взглянуть на содержимое GC вдоль хромосомы.

Вы извлекаете блоки 500000bp и вычисляете содержимое GC.

Вычислите, сколько блоков использовать.

numNT = numel(chr1.Data);
blockSize = 500000;
numBlocks = floor(numNT/blockSize);

Одним из способов помочь MATLAB эффективности при работе с большими наборами данных является «предварительное выделение» пространства для данных. Это позволяет MATLAB выделять достаточно пространство для всех данных, а не увеличивать массив маленькими фрагментами. Это ускорит процесс, а также защитит вас от проблем, связанных с тем, что данные становятся слишком большими для хранения. Для получения дополнительной информации о предварительном выделении массивов смотрите: http://www.mathworks.com/support/solutions/data/1-18150.html?solution=1-18150

Легкий способ предварительно выделить массив - использовать zeros функция.

ratio = zeros(numBlocks+1,1);

Закольцовывайте данные, ищущие C или G, а затем разделите это число на общее число A, T, C и G. Для запуска потребуется около минуты.

A = nt2int('A'); C = nt2int('C'); G = nt2int('G'); T = nt2int('T');

for count = 1:numBlocks
    % calculate the indices for the block
    start = 1 + blockSize*(count-1);
    stop = blockSize*count;
    % extract the block
    block = chr1.Data(start:stop);
    % find the GC and AT content
    gc = (sum(block == G | block == C));
    at = (sum(block == A | block == T));
    % calculate the ratio of GC to the total known nucleotides
    ratio(count) = gc/(gc+at);
end

Конечный блок меньше, поэтому относитесь к этому как к особому случаю.

block = chr1.Data(stop+1:end);
gc = (sum(block == G | block == C));
at = (sum(block == A | block == T));
ratio(end) = gc/(gc+at);

График содержимого GC для хромосомы Homo Sapiens 1

xAxis = [1:blockSize:numBlocks*blockSize, numNT];
plot(xAxis,ratio)
xlabel('Base pairs');
ylabel('Relative GC content');
title('Relative GC content of Homo Sapiens Chromosome 1')

Область в центре графика вокруг 140Mbp - это большая область Ns.

seqdisp(chr1.Data(140000000:140001000))

ans =

  17x71 char array

    '   1  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    '  61  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 121  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 181  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 241  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 301  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 361  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 421  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 481  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 541  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 601  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 661  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 721  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 781  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 841  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 901  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN'
    ' 961  NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN NNNNNNNNNN N                    '

Нахождение областей высокого газового Содержимого

Можно использовать find идентифицировать области с высоким содержимым GC.

indices = find(ratio>0.5);
ranges = [(1 + blockSize*(indices-1)), blockSize*indices];
fprintf('Region %d:%d has GC content %f\n',[ranges ,ratio(indices)]')

Region 500001:1000000 has GC content 0.501412
Region 1000001:1500000 has GC content 0.598332
Region 1500001:2000000 has GC content 0.539498
Region 2000001:2500000 has GC content 0.594508
Region 2500001:3000000 has GC content 0.568620
Region 3000001:3500000 has GC content 0.584572
Region 3500001:4000000 has GC content 0.548137
Region 6000001:6500000 has GC content 0.545072
Region 9000001:9500000 has GC content 0.506692
Region 10500001:11000000 has GC content 0.511386
Region 11500001:12000000 has GC content 0.519874
Region 16000001:16500000 has GC content 0.513082
Region 17500001:18000000 has GC content 0.513392
Region 21500001:22000000 has GC content 0.505598
Region 156000001:156500000 has GC content 0.504446
Region 156500001:157000000 has GC content 0.504090
Region 201000001:201500000 has GC content 0.502976
Region 228000001:228500000 has GC content 0.511946

Если вы хотите удалить временный файл, сначала необходимо очистить memmapfile объект.

clear chr1
delete(mmFilename)