In [16]:
%%bash

#Create an installation directory
mkdir src
cd src

#Take a clone of BioPerl from GitHub and change to the correct version:
git clone https://github.com/bioperl/bioperl-live.git
cd bioperl-live
git checkout bioperl-release-1-6-1
cd -

#Install the Ensembl git tools to allow easy management of the Ensembl repos.
git clone https://github.com/Ensembl/ensembl-git-tools.git
export PATH=$PWD/ensembl-git-tools/bin:$PATH

#Install the APIs that you need. You can install all the APIs using the git ensembl command:
git ensembl --clone api
/home/dhimmels/Documents/serg/rephetio/construct/gwas-catalog/ensembl/src
* Processing 'api'

* Working with module 'ensembl'
* Cloning from remote 'https://github.com/Ensembl/ensembl.git'
* Enabling git hooks

* Working with module 'ensembl-compara'
* Cloning from remote 'https://github.com/Ensembl/ensembl-compara.git'
* Enabling git hooks

* Working with module 'ensembl-funcgen'
* Cloning from remote 'https://github.com/Ensembl/ensembl-funcgen.git'
* Enabling git hooks

* Working with module 'ensembl-io'
* Cloning from remote 'https://github.com/Ensembl/ensembl-io.git'
* Enabling git hooks

* Working with module 'ensembl-variation'
* Cloning from remote 'https://github.com/Ensembl/ensembl-variation.git'
* Enabling git hooks

Cloning into 'bioperl-live'...
Note: checking out 'bioperl-release-1-6-1'.

You are in 'detached HEAD' state. You can look around, make experimental
changes and commit them, and you can discard any commits you make in this
state without impacting any branches by performing another checkout.

If you want to create a new branch to retain commits you create, you may
do so (now or later) by using -b with the checkout command again. Example:

  git checkout -b new_branch_name

HEAD is now at 83dcde3... * bump point version * sync minor doc changes
Cloning into 'ensembl-git-tools'...
Cloning into 'ensembl'...
Cloning into 'ensembl-compara'...
Cloning into 'ensembl-funcgen'...
Cloning into 'ensembl-io'...
Cloning into 'ensembl-variation'...
In [19]:
snps = ['rs1333049', 'rs2248462', 'rs225190', 'rs225212', 'rs2252865', 'rs2252996']
with open('SNP_list_sorted.txt', 'w') as write_file:
    write_file.write('\n'.join(sorted(snps)))
In [20]:
%%perl

use lib "src/bioperl-live";
use lib "src/ensembl/modules";
use lib "src/ensembl-compara/modules";
use lib "src/ensembl-variation/modules";
use lib "src/ensembl-funcgen/modules";

use strict;
use warnings;
use Bio::EnsEMBL::Registry;

open (IN, "<SNP_list_sorted.txt");
open (OUT, ">SNPs_in_LD.tsv");

my $reg = 'Bio::EnsEMBL::Registry';
$reg->load_registry_from_db(-host => 'ensembldb.ensembl.org',-user => 'anonymous');
my $va = $reg->get_adaptor( 'human', 'variation', 'variation' );

while (<IN>) {
    print $_;
    my $v = $va->fetch_by_name($_);
    if (!$v) {
        print 'unfetched ' . $_;
        next;
    }
    my @vfs = @{ $v->get_all_VariationFeatures };
    
    foreach my $vf (@vfs){        
        my $ld = $vf->get_all_LD_values;
        my @pops = @{ $vf->get_all_LD_Populations };
        my @ldvs = @{ $ld->get_variations };
        
        foreach my $pop (@pops) {
            
            if ($pop->name =~ /1000GENOMES/) {
        
                foreach my $ldv (@ldvs) {            
                    if ($ldv->stable_id ne $_) {
                        my @ldvfs = @{ $ldv->get_all_VariationFeatures };
                
                        foreach my $ldvf (@ldvfs) {
                            my @tvs = @{ $ldvf->get_all_TranscriptVariations };
                            my $r2 = $ld->get_r_square($vf, $ldvf, $pop);                                        
                    
                            foreach my $tv (@tvs) {
                                my $gene = $tv->transcript->get_Gene;
                            
                                if ($r2) {            
                                    print OUT $v->stable_id, "\t", $ldv->stable_id, "\t", $gene->external_name, "\t", $r2, "\t", $pop->name, "\n";
                                    }
                                }
                            }
                        }
                    }
                }
            }
        }
    }
rs1333049
unfetched rs1333049
rs2248462
unfetched rs2248462
rs225190
unfetched rs225190
rs225212
unfetched rs225212
rs2252865
unfetched rs2252865
rs2252996
UNIVERSAL->import is deprecated and will be removed in a future perl at src/bioperl-live/Bio/Tree/TreeFunctionsI.pm line 94.

-------------------- WARNING ----------------------
MSG: Binary file calc_genotypes not found. Please, read the ensembl-variation/C_code/README.txt file if you want to use LD calculation

FILE: Variation/DBSQL/LDFeatureContainerAdaptor.pm LINE: 674
CALLED BY: Variation/DBSQL/LDFeatureContainerAdaptor.pm  LINE: 607
Date (localtime)    = Mon Jun 22 18:03:46 2015
Ensembl API version = 80
---------------------------------------------------
In [ ]: