@@ -321,4 +321,51 @@ calculateBlastHomology <- function(gr, ref, db, anchorLength=150) {
321321 cl $ minOverlap <- pmin(cl $ leftOverlap , cl $ rightOverlap )
322322 return (cl )
323323}
324+ # ' Converts to breakend notation
325+ .toVcfBreakendNotationAlt = function (gr , insSeq = gr $ insSeq , ref = gr $ REF ) {
326+ assert_that(all(width(gr ) == 1 ))
327+ assert_that(! is.null(insSeq ))
328+ assert_that(all(insSeq != " " ))
329+ assert_that(! is.null(gr $ partner ))
330+ isBreakpoint = ! is.na(gr $ partner )
331+ breakendAlt = ifelse(as.character(strand(gr )) == " +" , paste0(gr $ insSeq , " ." ), paste0(" ." , gr $ insSeq ))
332+ gr $ partner [isBreakpoint ] = names(gr )[isBreakpoint ] # self partner to prevent errors
333+ partnergr = gr [gr $ partner ]
334+ partnerDirectionChar = ifelse(strand(partnergr ) == " +" , " ]" , " [" )
335+ breakpointAlt = ifelse(as.character(strand(gr )) == " +" ,
336+ paste0(ref , insSeq , partnerDirectionChar , seqnames(partnergr ), " :" , start(partnergr ), partnerDirectionChar ),
337+ paste0(partnerDirectionChar , seqnames(partnergr ), " :" , start(partnergr ), partnerDirectionChar , insSeq , ref ))
338+ return (ifelse(isBreakpoint , breakpointAlt , breakendAlt ))
339+ }
340+ # ' Converts the given breakpoint GRanges object to VCF format in breakend
341+ # ' notation.
342+ # '
343+ # ' @param gr breakpoint GRanges object. Can contain both breakpoint and single breakend SV records
344+ # '
345+ # '@export
346+ breakpointGRangesToVCF <- function (gr ) {
347+ if (is.null(gr $ insSeq )) {
348+ gr $ insSeq = rep(" " , length(gr ))
349+ }
350+ nominalgr = GRanges(seqnames = seqnames(gr ), ranges = IRanges(start = (end(gr ) + start(gr )) / 2 , width = 1 ))
351+ if (is.null(gr $ REF )) {
352+ gr $ REF = rep(" N" , length(gr ))
353+ }
354+ gr $ ALT [is.na(gr $ ALT )] = " "
355+ if (is.null(gr $ ALT )) {
356+ gr $ ALT = rep(" " , length(gr ))
357+ }
358+ gr $ ALT [is.na(gr $ ALT )] = " "
359+ gr $ ALT [gr $ ALT == " " ] = .toVcfBreakendNotationAlt(gr )[gr $ ALT == " " ]
360+ ciposstart = start(gr ) - start(nominalgr )
361+ ciposend = end(gr ) - end(nominalgr )
362+ vcf = VCF(rowRanges = nominalgr , collapsed = FALSE )
363+ fixeddf = data.frame (
364+ ALT = gr $ ALT ,
365+ REF = gr $ REF ,
366+ QUAL = gr $ QUAL ,
367+ FILTER = gr $ FILTER )
368+
369+ VCF(rowRanges = GRanges(), colData = DataFrame(), exptData = list (header = VCFHeader()), fixed = DataFrame(), info = DataFrame(), geno = SimpleList(), ... , collapsed = FALSE , verbose = FALSE
370+ }
324371
0 commit comments