Plotting relatedness2 output in R
0
0
Entering edit mode
11 months ago
Whirlingdaf ▴ 40

OK, I have a question that I feel has a quite simple solution, but my brain has just hit a wall...

I would like to plot the output from :

vcftools -relatedness2

into a plot where the samples are nodes and the edges show the values of relatedness. I would also like to be able to filter this relatedness to only show relatedness edges with values >0.01.

Example relatedness output looks something like:

INDV1, INDV2, RELATEDNESS_PHI

1,  4, 0.05096620

1,  2, 0.01837480

1,  3, -0.03895280

1,  6, 0.00689972
...

Where RELATEDNESS_PHI is the relatedness value. An example output would look something like:

enter image description here

Thanks for any help working past this mental block!

vcftools relatedness2 nodes R • 621 views
ADD COMMENT
0
Entering edit mode

And where exactly is now your mental block? The data you are showing corresponds to the default edge notation useable by e.g. Cytoscape or the geomnet package in R. For traversals through edge graphs (can even be bi-weighted) use dogdr.

ADD REPLY
0
Entering edit mode

Sorry just seeing this, I did end up trying to work through these data with igraph, so it is likely within my scripting in igraph itself... I had:

# Get the relatedness values and corresponding individual IDs

relatedness_values <- relatedness_data2$RELATEDNESS_PHI 
individuals1 <- relatedness_data$INDV1 
individuals2 <- relatedness_data$INDV2

# Combine all individual IDs
all_individuals <- unique(c(individuals1, individuals2))

# Create the graph 
relatedness_graph <- graph(edges = cbind(individuals1, individuals2), directed = FALSE)

# Set the vertex names to individual IDs
V(relatedness_graph)$name <- all_individuals

# Add relatedness values as edge attribute
E(relatedness_graph)$relatedness <- relatedness_values

# Subset edges greater than or equal to 0.08 and less than 0.5 
E(relatedness_graph)[[relatedness >= 0.08 & <= 0.5]] 

# Plot the graph
plot(relatedness_graph, vertex.label.color = "black", layout = layout_with_fr(g1))
ADD REPLY

Login before adding your answer.

Traffic: 2869 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6