<commit_before>#!/usr/bin/perl

=head1
load_genotypes_vcf_cxgn_postgres.pl - loading genotypes into cxgn databases, based on the load_cassava_snps.pl script by Naama.

=head1 SYNOPSIS
    perl bin/load_genotypes_vcf_cxgn_postgres.pl -H localhost -D fixture -U postgres -c VCF -o /tmp/transposevcf.txt -i /home/vagrant/Documents/cassava_subset_108KSNP_10acc.vcf -r /archive_path/ -g "test_pop_01" -p "test_project_01" -d "Diversity study" -y 2016 -l "test_location" -n "IGD" -b "accession" -m "test_protocol_01_new" -k "protocol description" -q "Manihot esculenta" -e "IITA" -s -u nmorales -f "Mesculenta_511_v7" -B /tmp/SQLCOPY.csv -A

    If you are loading a "transposed VCF" use -c transposedVCF otherwise for a normal VCF use -c VCF. When using a normal VCF, give a temporary file using -o for where this script will transpose your VCF. VCF are transposed for speed.
    To use an existing project (not create a new project name entry), use -h project_id
    To use an existing protocol (not create a new nd_protocol name entry), use -j protocol_id

=head1 COMMAND-LINE OPTIONS
  ARGUMENTS
 -H host name (required) e.g. "localhost"
 -D database name (required) e.g. "cxgn_cassava"
 -U database username (required)
 -c VCF file type. transposedVCF or VCF
 -o temporary file for transposing a VCF. whenever a VCF is used, it is transposed for speed (unless option w is flagged).
 -i path to infile (required)
 -r archive path (required)
 -u username in database (required)
 -p project name (required) e.g. "SNP genotyping 2012 Cornell Biotech".  Will be found or created in Project table.
 -y project year (required) e.g. "2012".  Will be saved as a Projectprop.
 -d project description (required) e.g. "Diversity study"
 -n genotype facility name (required) e.g. "igd"
 -g population name (required) e.g. "NaCRRI training population"
 -b observation unit type name (required) e.g. "tissue_sample" or "accession"
 -e breeding program name (required) e.g. "IITA"
 -s include igd numbers in sample names
 -m protocol name (required) e.g. "GBS ApeKI Cassava genome v6"
 -k protocol description (required)
 -l location name (required) e.g. "Cornell Biotech".  Will be found or created in NdGeolocation table.
 -q organism species name (required) e.g. "Manihot esculenta".
 -f reference genome name (required) e.g. "Mesculenta_511_v7"
 -B temporary file where the SQL COPY file is written.

 -h project_id (Will associate genotype data to an existing project_id)
 -j protocol_id (Will associate genotype data to an existing nd_protocol_id)

  FLAGS
 -x delete old genotypes for accessions that have new genotypes
 -a add accessions that are not in the database
 -z if accession names include an IGD number. Accession names are in format 'accession_name:IGD_number'. The IGD number will be parsed and stored as a genotypeprop.
 -t Test run . Rolling back at the end.
 -w in the case that you have uploaded a normal VCF and you do not want to transpose it (because the transposition is memory intensive), use this flag
 -A accept warnings and continue with the storing. warnings are whether the samples already have genotype scores for a specific protocol/project

=head1 DESCRIPTION
This script loads genotype data into the Chado genotype table it encodes the genotype + marker name in a json format in the genotyope.uniquename field for easy parsing by a Perl program. The genotypes are linked to the relevant stock using nd_experiment_genotype. Each column in the spreadsheet, which represents a single accession (stock) is stored as a single genotype entry and linked to the stock via nd_experiment_genotype. Stock names are stored in the stock table if cannot be found, and linked to a population stock with the name supplied in opt_g. Map details (chromosome, position, ref, alt, qual, filter, info, and format) are stored in json format in the protocolprop table.

This script mimics exactly the "online process" in SGN::Controller::AJAX::GenotypesVCFUpload->upload_genotype_verify

Can use a transposedVCF or normal VCF

=head1 AUTHOR
 <NAME> (<EMAIL>)
 <NAME> <<EMAIL>>
=cut

use strict;

use Getopt::Std;
use Data::Dumper;
use JSON::Any;
use JSON::PP;
use Carp qw /croak/ ;
use Try::Tiny;
use Pod::Usage;

use Bio::Chado::Schema;
use CXGN::Metadata::Schema;
use CXGN::Phenome::Schema;
use CXGN::People::Person;
use CXGN::DB::InsertDBH;
use Sort::Versions;
use SGN::Model::Cvterm;
use CXGN::Genotype::StoreVCFGenotypes;
use DateTime;
use CXGN::UploadFile;
use File::Basename qw | basename dirname|;
use CXGN::Genotype::Protocol;
use CXGN::Genotype::ParseUpload;

our ($opt_H, $opt_D, $opt_U, $opt_c, $opt_o, $opt_r, $opt_i, $opt_t, $opt_p, $opf_f, $opt_y, $opt_g, $opt_a, $opt_x, $opt_v, $opt_s, $opt_m, $opt_k, $opt_l, $opt_q, $opt_z, $opt_u, $opt_b, $opt_n, $opt_s, $opt_e, $opt_f, $opt_d, $opt_h, $opt_j, $opt_w, $opt_A, $opt_B);

getopts('H:U:i:r:u:c:o:tD:p:y:g:axsm:k:l:q:zf:d:b:n:se:h:j:wAB:');

if (!$opt_H || !$opt_U || !$opt_D || !$opt_c || !$opt_i || !$opt_p || !$opt_y || !$opt_m || !$opt_k || !$opt_l || !$opt_q || !$opt_r || !$opt_u || !$opt_f || !$opt_d || !$opt_b || !$opt_n || !$opt_e || !$opt_B) {
    pod2usage(-verbose => 2, -message => "Must provide options -H (hostname), -D (database name), -U (database username), -c VCF file type (transposedVCF or VCF), -i (input file), -r (archive path), -p (project name), -y (project year), -l (location name of project), -m (protocol name), -k (protocol description), -q (organism species), -u (database username), -f (reference genome name), -d (project description), -b (observation unit type name), -n (genotype facility name), -e (breeding program name), -B (temp file where SQL COPY is written)\n");
}

if ($opt_c ne 'transposedVCF' && $opt_c ne 'VCF') {
    die "Not a valid option c\n";
}

if ($opt_c eq 'VCF '&& !$opt_o) {
    die "When uploading a VCF e.g. option c is VCF, you must give a temporary file using option o, so that this script can transpose your file before loading. All VCF are transposed for speed of loading.\n";
}

my $file = $opt_i;
my $protocol_id = $opt_j;
my $organism_species = $opt_q;
my $reference_genome_name = $opt_f;
my $obs_type = $opt_b;
my $add_accessions = 0;
if ($opt_a){
    $add_accessions = 1;
}
my $include_igd_numbers = 0;
if ($opt_z){<commit_msg>opt s and opt z the same. left only opt z<commit_after>#!/usr/bin/perl

=head1
load_genotypes_vcf_cxgn_postgres.pl - loading genotypes into cxgn databases, based on the load_cassava_snps.pl script by Naama.

=head1 SYNOPSIS
    perl bin/load_genotypes_vcf_cxgn_postgres.pl -H localhost -D fixture -U postgres -c VCF -o /tmp/transposevcf.txt -i /home/vagrant/Documents/cassava_subset_108KSNP_10acc.vcf -r /archive_path/ -g "test_pop_01" -p "test_project_01" -d "Diversity study" -y 2016 -l "test_location" -n "IGD" -b "accession" -m "test_protocol_01_new" -k "protocol description" -q "Manihot esculenta" -e "IITA" -z -u nmorales -f "Mesculenta_511_v7" -B /tmp/SQLCOPY.csv -A

    If you are loading a "transposed VCF" use -c transposedVCF otherwise for a normal VCF use -c VCF. When using a normal VCF, give a temporary file using -o for where this script will transpose your VCF. VCF are transposed for speed.
    To use an existing project (not create a new project name entry), use -h project_id
    To use an existing protocol (not create a new nd_protocol name entry), use -j protocol_id

=head1 COMMAND-LINE OPTIONS
  ARGUMENTS
 -H host name (required) e.g. "localhost"
 -D database name (required) e.g. "cxgn_cassava"
 -U database username (required)
 -c VCF file type. transposedVCF or VCF
 -o temporary file for transposing a VCF. whenever a VCF is used, it is transposed for speed (unless option w is flagged).
 -i path to infile (required)
 -r archive path (required)
 -u username in database (required)
 -p project name (required) e.g. "SNP genotyping 2012 Cornell Biotech".  Will be found or created in Project table.
 -y project year (required) e.g. "2012".  Will be saved as a Projectprop.
 -d project description (required) e.g. "Diversity study"
 -n genotype facility name (required) e.g. "igd"
 -g population name (required) e.g. "NaCRRI training population"
 -b observation unit type name (required) e.g. "tissue_sample" or "accession"
 -e breeding program name (required) e.g. "IITA"
 -m protocol name (required) e.g. "GBS ApeKI Cassava genome v6"
 -k protocol description (required)
 -l location name (required) e.g. "Cornell Biotech".  Will be found or created in NdGeolocation table.
 -q organism species name (required) e.g. "Manihot esculenta".
 -f reference genome name (required) e.g. "Mesculenta_511_v7"
 -B temporary file where the SQL COPY file is written.

 -h project_id (Will associate genotype data to an existing project_id)
 -j protocol_id (Will associate genotype data to an existing nd_protocol_id)

  FLAGS
 -x delete old genotypes for accessions that have new genotypes
 -a add accessions that are not in the database
 -z if accession names include an IGD number. Accession names are in format 'accession_name:IGD_number'. The IGD number will be parsed and stored as a genotypeprop.
 -t Test run . Rolling back at the end.
 -w in the case that you have uploaded a normal VCF and you do not want to transpose it (because the transposition is memory intensive), use this flag
 -A accept warnings and continue with the storing. warnings are whether the samples already have genotype scores for a specific protocol/project

=head1 DESCRIPTION
This script loads genotype data into the Chado genotype table it encodes the genotype + marker name in a json format in the genotyope.uniquename field for easy parsing by a Perl program. The genotypes are linked to the relevant stock using nd_experiment_genotype. Each column in the spreadsheet, which represents a single accession (stock) is stored as a single genotype entry and linked to the stock via nd_experiment_genotype. Stock names are stored in the stock table if cannot be found, and linked to a population stock with the name supplied in opt_g. Map details (chromosome, position, ref, alt, qual, filter, info, and format) are stored in json format in the protocolprop table.

This script mimics exactly the "online process" in SGN::Controller::AJAX::GenotypesVCFUpload->upload_genotype_verify

Can use a transposedVCF or normal VCF

=head1 AUTHOR
 <NAME> (<EMAIL>)
 <NAME> <<EMAIL>>
=cut

use strict;

use Getopt::Std;
use Data::Dumper;
use JSON::Any;
use JSON::PP;
use Carp qw /croak/ ;
use Try::Tiny;
use Pod::Usage;

use Bio::Chado::Schema;
use CXGN::Metadata::Schema;
use CXGN::Phenome::Schema;
use CXGN::People::Person;
use CXGN::DB::InsertDBH;
use Sort::Versions;
use SGN::Model::Cvterm;
use CXGN::Genotype::StoreVCFGenotypes;
use DateTime;
use CXGN::UploadFile;
use File::Basename qw | basename dirname|;
use CXGN::Genotype::Protocol;
use CXGN::Genotype::ParseUpload;

our ($opt_H, $opt_D, $opt_U, $opt_c, $opt_o, $opt_r, $opt_i, $opt_t, $opt_p, $opf_f, $opt_y, $opt_g, $opt_a, $opt_x, $opt_v, $opt_m, $opt_k, $opt_l, $opt_q, $opt_z, $opt_u, $opt_b, $opt_n, $opt_e, $opt_f, $opt_d, $opt_h, $opt_j, $opt_w, $opt_A, $opt_B);

getopts('H:U:i:r:u:c:o:tD:p:y:g:axsm:k:l:q:zf:d:b:n:e:h:j:wAB:');

if (!$opt_H || !$opt_U || !$opt_D || !$opt_c || !$opt_i || !$opt_p || !$opt_y || !$opt_m || !$opt_k || !$opt_l || !$opt_q || !$opt_r || !$opt_u || !$opt_f || !$opt_d || !$opt_b || !$opt_n || !$opt_e || !$opt_B) {
    pod2usage(-verbose => 2, -message => "Must provide options -H (hostname), -D (database name), -U (database username), -c VCF file type (transposedVCF or VCF), -i (input file), -r (archive path), -p (project name), -y (project year), -l (location name of project), -m (protocol name), -k (protocol description), -q (organism species), -u (database username), -f (reference genome name), -d (project description), -b (observation unit type name), -n (genotype facility name), -e (breeding program name), -B (temp file where SQL COPY is written)\n");
}

if ($opt_c ne 'transposedVCF' && $opt_c ne 'VCF') {
    die "Not a valid option c\n";
}

if ($opt_c eq 'VCF '&& !$opt_o) {
    die "When uploading a VCF e.g. option c is VCF, you must give a temporary file using option o, so that this script can transpose your file before loading. All VCF are transposed for speed of loading.\n";
}

my $file = $opt_i;
my $protocol_id = $opt_j;
my $organism_species = $opt_q;
my $reference_genome_name = $opt_f;
my $obs_type = $opt_b;
my $add_accessions = 0;
if ($opt_a){
    $add_accessions = 1;
}
my $include_igd_numbers = 0;
if ($opt_z){
{'fraction_non_alphanumeric': 0.11272115923464894, 'fraction_numerical': 0.006724841886158034, 'mean_word_length': 4.447884866986481, 'pattern_counts': {'":': 0, '<': 13, '<?xml version=': 0, '>': 19, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '4110131', 'n_tokens_mistral': 4430, 'n_tokens_neox': 4215, 'n_words': 1702}