Если .. еще, если .. еще цикл в R - PullRequest
Новая
       11

Если .. еще, если .. еще цикл в R

0 голосов
/ 18 октября 2018

У меня довольно длинный цикл в R, который я пытаюсь использовать для выравнивания двух фреймов данных (dosage и ldpred) по значениям, которые они имеют в определенных столбцах.Я заранее извинюсь за стену текста, которая будет следовать, но это выглядит как довольно сложная вещь (возможно, с очень простым решением).

Цикл должен принимать значение строки snp из ldpred и найдите ту же строку в dosage, а затем укажите конкретную строку из dosage.Затем предполагается использовать значения двух других столбцов в ldpred и сравнить их с соответствующими столбцами в dosage.Если значения совпадают, то в новом столбце должно быть 1.Если они совпадают, но находятся в противоположных столбцах, он должен вывести -1.И это нормально.

Однако, сложность заключается в том, что он также должен затем переключать значения ldpred на другие значения и затем повторять ту же проверку, что и выше.

Если все это не возвращает 1 или -1 (т. Е. Значения просто не совпадают по какой-либо причине), то предполагается вернуть 0.

РЕДАКТИРОВАТЬ по запросу, содержащему примеры моих данныхи вывод:

Дозировка: 

chr snp a1 a2 p-value
1 rs1234 A G 0.05
2 rs2345 C T 0.03
3 rs5555 G T 0.001
4 rs9876 C G 0.02

LDpred: 

chr sid nt1 nt2 beta OUTPUT
1 rs1234 A G 0.001 1
2 rs2345 T C 0.002 -1
3 rs5555 C A 0.003 1
4 rs9876 CC GG 0.004 0

Надеюсь, это немного прояснит ситуацию.Я пытаюсь найти значение SNP в LDpred, найти соответствующее значение SNP в дозировке, а затем сравнить значения nt1 с a1 и nt2 с a2. / edit

Вот скрипт: 

for (line in 1:nrow(ldpred)){

  # Input rsID and genotype of specific line of LDpred file

  snp_ld = ldpred$sid[line]
  ref_ld = ldpred$nt1[line]
  alt_ld = ldpred$nt2[line]

  # Obtain opposing line from dosage file using rsID

  genotype = subset(dosage, snp == snp_ld)

  # Extract dosage file genotypes from dosage line

  ref_gen = genotype$a1
  alt_gen = genotype$a2

  if (ref_ld == ref_gen && alt_ld == alt_gen){

    # If alleles in both files match, return 1

    ldpred$matched[line] = 1

    }  else if (ref_ld == alt_gen && alt_ld == ref_gen){

      # If alleles in both files are exact opposites, return -1

      ldpred$matched[line] = -1

      }  else{

        # Make sure that files aren't using alternate strands
        # Switch alleles to opposing strand using switch_strand function

        ref_ld_switched = switch_strand(ref_ld)
        alt_ld_switched = switch_strand(alt_ld)

        if (ref_ld_switched == ref_gen && alt_ld_switched == alt_gen){

          # If new switched alleles match, return 1

          ldpred$matched[line] = 1

          }  else if (ref_ld == alt_gen && alt_ld == ref_gen){

            # If new switched alleles are opposites, return -1 

            ldpred$matched[line] = -1

          }  else {

            # If the alleles do not match then return 0 for QC

            ldpred$matched = 0

          }
      }
}

У меня изначально было довольно много проблем с использованием цикла for .. if .. else и соответствующего фигурногоскобки, но я думаю Я отсортировал это (однако, если кто-то обнаружит какие-либо ошибки, я буду признателен за информацию).Теперь, хотя я получаю сообщение об ошибке 

Error in if (ref_ld == ref_gen && alt_ld == alt_gen) { :   
  missing value where TRUE/FALSE needed

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

Любая помощь приветствуется!

Ответы [ 2 ]

0 голосов
/ 18 октября 2018

R очень хорошо работает с векторами и позволяет легко индексировать ваши данные, используя скобки.Лучше воспользоваться этим и по возможности избегать использования петель.Вот альтернативы вашему коду: 

dosage <- read.table(text =
  'chr snp a1 a2 p-value
  1 rs1234 A G 0.05
  2 rs2345 C T 0.03
  3 rs5555 G T 0.001
  4 rs9876 C G 0.02',
  header = TRUE, stringsAsFactors = FALSE)

ldpred <- read.table(text =
  'chr sid nt1 nt2 beta
  1 rs1234 A G 0.001
  2 rs2345 T C 0.002
  3 rs5555 C A 0.003',
  header = TRUE, stringsAsFactors = FALSE)
# I removed the last line to show what happens when you have different sizes

mergedData <- merge(ldpred, dosage, by.x = c('chr','sid'), by.y = c('chr','snp'), all = TRUE)

mergedData$OUTPUT <- 0
mergedData$OUTPUT[mergedData$nt1 == mergedData$a1 & mergedData$nt2 == mergedData$a2] <- 1
mergedData$OUTPUT[mergedData$nt1 == mergedData$a2 & mergedData$nt2 == mergedData$a1] <- -1
mergedData$OUTPUT[apply(mergedData[,c('nt1','nt2','a1','a2')], 1, anyNA)] <- NA

> mergedData
  chr    sid  nt1  nt2  beta a1 a2 p.value OUTPUT
1   1 rs1234    A    G 0.001  A  G   0.050      1
2   2 rs2345    T    C 0.002  C  T   0.030     -1
3   3 rs5555    C    A 0.003  G  T   0.001      0
4   4 rs9876 <NA> <NA>    NA  C  G   0.020     NA

merge () гарантирует, что ваши данные выровнены для выполнения векторного сравнения, вы можете сделать это другими способами.Вы также можете удалить опцию all = TRUE, чтобы у вас были только совпадающие строки, избегая NA.Обычно я этого не делаю, потому что АН - это информация.

0 голосов
/ 18 октября 2018

Я думаю, что задним числом это был неправильный форум, чтобы задавать такой конкретный вопрос, который по существу был связан с отладкой.

Для потомков: я определил определил ответ на свой вопрос было на самом деле, потому что файлы данных были разных размеров.Таким образом, в итоге переменная genotype оказалась нулевым значением, и когда она была введена в первый этап цикла if, она выдала ошибку «пропущенные значения».

...