library(readr)
library(dplyr)
# --- 1. Cargar los datos ---
file_path = "FAM8.xlsx - CSC_32553.varfile.ex.EDITED.csv"
df = read_csv(file_path)
names(df)## [1] "#CHR" "START" "END" "REF"
## [5] "ALT" "ZYG" "FILTER" "QUAL"
## [9] "RR" "VR" "TR" "RATIO"
## [13] "GQ" "PL" "VCF_INFO" "BAND"
## [17] "SNP17NF" "ANNOTATION" "GENE" "GENE_NAME"
## [21] "ACCESSION" "EXON" "NT_CHANGE" "AA_CHANGE"
## [25] "ALL_TRX" "TGP_FREQ" "ESP_FREQ" "GRANTHAM_SC"
## [29] "PHASTCONS" "GERP_ESP" "GERP" "PHYLOP_SC"
## [33] "PHYLOP_PRED" "LRT_SCORE" "LRT_PRED" "SIFT_SCORE"
## [37] "SIFT_PRED" "PPH2_SCORE" "PPH2_PRED" "MTT_SCORE"
## [41] "MTT_PRED" "ESP_PPH2" "CHIMP" "EVS_ALLELES"
## [45] "ESP_GENOTYPES" "EVE_ALT_FREQ" "EVE_HOM" "Proband:ZYG"
## [49] "Father:ZYG" "Mother:ZYG" "UASib:ZYG"## [1] 22823 51df <- df %>%
mutate(
TGP_FREQ = na_if(TGP_FREQ, "na"),
ESP_FREQ = na_if(ESP_FREQ, "na")
)
df <- df %>%
mutate(
TGP_FREQ = as.numeric(na_if(TGP_FREQ, "na")),
ESP_FREQ = as.numeric(na_if(ESP_FREQ, "na"))
)Identificacion de la variante
#CHR → Cromosoma en donde esta la variante (1-22,X,Y,MT)
START → Posición inicial en el cromosoma (coordenada genómica).
END → Posición final (para SNVs suele ser igual a START, para indels puede diferir).
REF → Alelo de referencia (según el genoma de referencia, ej. GRCh38).
ALT → Alelo alternativo observado en tu muestra.
ZYG → Cigocidad en la muestra del calling original (puede variar con las columnas por individuo).
Calidad del llamado
FILTER → Estado del filtro
QUAL → Puntaje de calidad de la llamada (Phred-scale, más alto = más confiable).
RR → Read count de alelos de referencia (cuántos reads apoyan al REF).
VR → Read count de alelos variantes (cuántos reads apoyan al ALT).
TR → Total de lecturas en esa posición (RR+VR).
RATIO → Proporción de alelo variante (VR/TR).
GQ → Genotype Quality (calidad del genotipo asignado, en escala Phred).
PL → Genotype likelihoods (valores de probabilidad para hom REF / het / hom ALT).
Informacion de VCF extendida
VCF_INFO → Campo extra con anotaciones del VCF
(DP, MQ, QD, etc.).
BAND → Región cromosómica (ej. 1q41).
SNP17NF → ID de la variante en una base de datos de SNPs.
Anotacion funcional
ANNOTATION → Tipo de cambio funcional (ej.
nonsynonymous SNV, synonymous,
stopgain, frameshift).
GENE → Símbolo del gen afectado (ej. MARCH1).
GENE_NAME → Nombre extendido del gen.
ACCESSION → Transcript ID en RefSeq (ej. NM_022746).
EXON → Exón afectado.
NT_CHANGE → Cambio a nivel de nucleótido (ej. c.A493G).
AA_CHANGE → Cambio a nivel de aminoácido (ej. p.T165A).
ALL_TRX → Anotación de todos los transcritos (no solo el principal).
Frecuencia poblacional
TGP_FREQ → Frecuencia alélica en 1000 Genomes Project.
ESP_FREQ → Frecuencia alélica en Exome Sequencing Project.
EVE_ALT_FREQ → Otras bases de datos de frecuencia de alelos (ExAC, gnomAD, etc. según el pipeline).
EVE_HOM → Número de homocigotos observados en esas bases.
EVS_ALLELES / ESP_GENOTYPES → Información de genotipos en cohortes de referencia.
Predictores funcionales
(miden el efecto de variantes missense o no sinonimas)
GRANTHAM_SC → Distancia química entre aminoácidos (alto valor = cambio más drástico).
PHASTCONS → Conservación filogenética (0 = no conservado, 1 = altamente conservado).
GERP / GERP_ESP → Puntuación de conservación evolutiva (alto = función importante).
PHYLOP_SC → Conservación evolutiva (positivo = conservado, negativo = acelerado).
PHYLOP_PRED → Predicción (conservado / neutro).
LRT_SCORE / LRT_PRED → Likelihood Ratio Test, predice si la variante es funcional o neutral.
SIFT_SCORE / SIFT_PRED → Predice tolerado o dañino basado en conservación de aminoácidos.
PPH2_SCORE / PPH2_PRED → PolyPhen-2, clasifica como benigno, posiblemente dañino o probablemente dañino.
MTT_SCORE / MTT_PRED → Otro predictor de patogenicidad (dependiendo del pipeline).
ESP_PPH2 → Predicciones integradas de PolyPhen-2 en ESP.
Evolutivo/Comparativo
CHIMP → Alelo en chimpancé (para inferir ancestralidad).
EVE_ALT_FREQ → Frecuencia alternativa en bases evolutivas.
Genotipos familiares
Proband:ZYG → Cigocidad de la variante en el
paciente (ej. het, hom).
Father:ZYG → Cigocidad en el padre.
Mother:ZYG → Cigocidad en la madre.
UASib:ZYG → Cigocidad en la hermana no afectada.
Se puede ver:
De novo: Variante presente en el probando, ausente en ambos padres.
Recesivo homocigoto: Proband hom ALT, padres heterocigotos.
Compuesto heterocigoto: Proband heterocigoto para dos variantes diferentes en el mismo gen, cada variante heredada de un progenitor distinto.
#filtrado por calidad de QC, variantes sinonimas, y MAF <
df_filt = df %>%
# Mantienes el filtro de calidad y exclusión de sinónimas
filter(QUAL > 30) %>%
filter(!ANNOTATION %in% c("synonymous SNV", "nonframeshift insertion",
"nonframeshift deletion", "ncRNA_splicing")) %>%
# Filtros de frecuencia alélica (MAF)
filter(is.na(TGP_FREQ) | TGP_FREQ <= 0.01) %>%
filter(is.na(ESP_FREQ) | ESP_FREQ <= 0.01) %>%
filter(is.na(EVE_ALT_FREQ) | EVE_ALT_FREQ <= 0.01) %>%
# Filtro de calidad de genotipo
filter(GQ >= 20 | is.na(GQ)) %>%
# Filtros para predictores funcionales
filter(is.na(SIFT_PRED) | SIFT_PRED %in% c("damaging")) %>%
filter(is.na(LRT_PRED) | LRT_PRED %in% c("deleterious")) %>%
filter(is.na(MTT_PRED) | MTT_PRED %in% c("disease_causing"))# De novo: presente en probando, ausente en ambos padres y hermana
denovo_variants <- df_filt %>%
filter(`Proband:ZYG` %in% c("Proband:het", "Proband:hom")) %>%
filter(`Father:ZYG` %in% c("Father:na")) %>%
filter(`Mother:ZYG` %in% c("Mother:na")) %>%
filter(`UASib:ZYG` %in% c("UASib:na"))
# Lista de genes candidatos de novo
denovo_genes <- unique(denovo_variants$GENE)
print(denovo_genes)## [1] "CPS1" "SLC9B1" "TUBB3"# Hom recesivos: probando hom ALT, padres heterocigotos
recessive_hom <- df_filt %>%
filter(`Proband:ZYG` == "Proband:hom") %>%
filter(`Father:ZYG` == "Father:het") %>%
filter(`Mother:ZYG` == "Mother:het")
recessive_hom_genes <- unique(recessive_hom$GENE)
print(recessive_hom_genes)## character(0)#Heterocigosidad compuesta
# Filtra las variantes donde el probando es heterocigoto
comphet_candidates = df_filt %>%
filter(`Proband:ZYG` == "Proband:het")
# Identifica genes donde el probando tiene dos variantes diferentes,
# una heredada de cada padre
comphet_genes = comphet_candidates %>%
# Agrupa por gen para analizar las variantes juntas
group_by(GENE) %>%
# Asegúrate de que hay al menos 2 variantes en el mismo gen
filter(n() >= 2) %>%
# Verifica que hay una variante con el patrón de herencia paterno (Father:het)
# y otra variante con el patrón de herencia materno (Mother:het)
filter(any(`Father:ZYG` == "Father:het" & `Mother:ZYG` == "Mother:na") &
any(`Mother:ZYG` == "Mother:het" & `Father:ZYG` == "Father:na")) %>%
# Selecciona el nombre del gen
pull(GENE) %>%
unique()
print(comphet_genes)## character(0)## [1] "CPS1" "SLC9B1" "TUBB3"TUBB3 es el mas coherente con el fenotipo observado y descrito para la familia 8.
Phenotype: Microcephaly; Cerebellar malformation; Seizures, febrile; Global developmental delay; Vision impairment, central; Hypotonia; Dystonia, generalized; Hyperreflexia; Weak cry
Ademas tiene sentido que sea una variante de novo dado a que ni los padres y hermana estan afectados con la enfermedad.
## [,1]
## #CHR "16"
## START "90001392"
## END "90001392"
## REF "C"
## ALT "T"
## ZYG "het"
## FILTER "VQSRTrancheSNP99.00to99.90"
## QUAL "1090.11"
## RR "53"
## VR "42"
## TR "95"
## RATIO "0.4421053"
## GQ "99"
## PL "113001536"
## VCF_INFO "AC=1;AF=0.500;AN=2;BaseQRankSum=0.393;ClippingRankSum=1.92;DP=95;FS=3.935;MQ=58.98;MQ0=0;MQRankSum=0.378;QD=11.47;ReadPosRankSum=0.738;VQSLOD=-4.041e-01;culprit=DP;;"
## BAND "16q24.3"
## SNP17NF "na"
## ANNOTATION "nonsynonymous SNV"
## GENE "TUBB3"
## GENE_NAME "tubulin, beta 3 class III"
## ACCESSION "NM_006086"
## EXON "exon4"
## NT_CHANGE "c.C533T"
## AA_CHANGE "p.T178M"
## ALL_TRX "TUBB3:NM_006086:exon4:c.C533T:p.T178M,TUBB3:NM_001197181:exon4:c.C317T:p.T106M,"
## TGP_FREQ NA
## ESP_FREQ NA
## GRANTHAM_SC "na"
## PHASTCONS "na"
## GERP_ESP "na"
## GERP "4.52"
## PHYLOP_SC "0.997528"
## PHYLOP_PRED "conserved"
## LRT_SCORE "0.999983"
## LRT_PRED "deleterious"
## SIFT_SCORE "0.917204"
## SIFT_PRED NA
## PPH2_SCORE "0.8"
## PPH2_PRED "possibly_damaging"
## MTT_SCORE "0.999782"
## MTT_PRED "disease_causing"
## ESP_PPH2 "na"
## CHIMP "na"
## EVS_ALLELES "na"
## ESP_GENOTYPES "0"
## EVE_ALT_FREQ "0"
## EVE_HOM "0"
## Proband:ZYG "Proband:het"
## Father:ZYG "Father:na"
## Mother:ZYG "Mother:na"
## UASib:ZYG "UASib:na"