cannin/biostars_qa huggingface dataset:
Question: I wish to process orthomcl output clusters with perl.

I have general idea about how the output looks like. but have no Linux resources right now.

Can anybody kindly run orthomcl on my data. (I will keep the input as small as possible just for testing purpose) and deliver the cluster file?

Default parameters (identity, coverage) will do.

Even results for this **small** input will suffice:

```
>1|307346730|Prevotella_amnii_CRIS_21A-A|31#
MLCFANDAIKEGYHVLGGDNYLMLDGNFVLEPKRINKSSGTSHIVISFIHKSNLYAGLYHIQMDLK
>1|306346738|Prevotella_amnii_CRIS_21A-A|35#
MLCFANDAIKEGYHVLGGDNYLMLDGNFVLEPKRINKSSGTSHIVISFIHKSNLYAGLYHIQMDLK
>1|307346705|Prevotella_amnii_CRIS_21A-A|32#
MKALKKIFNWYFSANALPYWVILIIDGIICYLSGLIVFWLSCRGNFISLTNFKSISCSIIIYMIFVLIGIRAFKTYSGIVRYSSFVDIRRVFFAILLSLIIAEAMHYVIYNWRFAGIEFYPMNGKQIFVMYI
>1|307346704|Prevotella_amnii_CRIS_21A-A|32#
MFSGIIEEMATIVAIKHEKENIHLTLTCSFVNELGIDQSISHNGVCLTVVCIEGERYTVTAMRETIERSNLGLLKVGDSVNVERSMIMNGRLDGHIVQGHVDETATCVDIKDANGSTYYTFEYPFNKAMAARGYFTVDKGSITVNGVSLTVCQPTDNTFTVAIIPYTKENTNFKYIGVGTIVNIEFDILGKYIARLQSIIK
>1|307346728|Prevotella_amnii_CRIS_21A-A|34#
MKQKYILLTVLGVALFSLTGCKKEKHTHDIITKIAPKPKLPKGPIELQGFTYTKKVDWLGTQYIITIIRKADKGLPLAYDEDGKKYYDNKIIINIKREDGTKFYQKEFTKKDFIQFSDEQYGKNGALLGFMFDSVYLKTLRFGASIGSPAPSSDEYIPIDIVIDNLGHLSMTKAQTLDTDGDQQKDKNNSKATEKAKDTEEEGV
>2|270333715|Prevotella_bergensis_DSM_17361|47#
MEIIAVESQAYQELIDRLNRIEQYVERTSRLIQDIDDELEMTTKDLIGTLNVSESTLYRWRKKQLVRYRYTEGGDVRYFFKSIVIATKCNRLRVSGMRNDEILGRLNRFKDNLIMSSCINPKNRQL
>2|270333716|Prevotella_bergensis_DSM_17361|45#
MIEKEQILLLTQGGLNVFSHFLGFEVNLHRNFRSPFYDDRRASCHIYYDRKTSSYKFYDHGDTTYSGDCFWFVATLRNLNLKTSFPEVLETIVQELGLYSLCDGEKHSSHITSAYKKTIVPTPKADITKWTEERPYSFEIQPFDDGLLNYWGHYGIHEDTLRRFRVRSLKRYESVSAEGKKFELYSSPTEPMFAYIGNGYVKIYRPHSPKIRFLYGGRMPATYCFGMEQIPAKGDMLFITGGEKDVLSLYAHGFNAICFNSETA
>2|270333717|Prevotella_bergensis_DSM_17361|42#
MPYIDYKTEDTWQKRLFDKLISVEDKLDRLLILQEQSVDTTVHPPLKPEYLDIIDVSKILKVEQKTIYNWVWAGKIPYLKANGRLLFLREEIDKMVRKRDDW
>2|270333718|Prevotella_bergensis_DSM_17361|38#
MANKESRSIQEQIDLLKHRGMIIEDEEFAHLHLSHISYYRLKGYWWDMQTDKERHIFKNDANFKDVIARYFFDKELRLILFDAIEAIEIALRTKMIYHLSQSYGGLYYMDKGLFNNEELQQQHIHDLMGEFM
```

Answer: Hi,

If you are interested, you can run it yourself on you local machine. [Here][1] is the virtual box image with the pre-installed orthomcl. If you are in windows, you can get it started in less than 5 minutes. There is an accompanying guide as well. I am the author of this guide/VM, so if you have any questions, I will be glad to answer!

[1]: https://wiki.itap.purdue.edu/display/BioCoreFacility/Phylogenomics+workshop
Question: Hi, I was mirroring the UCSC genomebrowser in Ubuntu 12.04 following the steps mentioned [here][1]. The very last step was to download the "gbdb" data using the command:

rsync -avzP --delete --max-delete=20 rsync://hgdownload.cse.ucsc.edu/gbdb/mm9/ /home/genomebrowser/gbdb/mm9/

It ran for a long time, and after downloading 139GB in my `/home` partition, my OS complains since that overflows that partition size (that I misjudged at first by reading the tutorial where the author states only 62GB!). So can you tell me what is the size of the download from `ucsc/gbdb/mm9`? Or how can I get that information about that rsync server? Thanks.

[1]: http://enotacoes.wordpress.com/2009/09/03/installing-a-minimal-ucsc-mirror-in-ubuntu-jaunty-64-bits/

Answer: I just got the answer by using the "rsync" command:

rsync --size-only --recursive rsync://hgdownload.cse.ucsc.edu/gbdb/mm9/

And the answer is, currently (today 1 July, 2012, 1:50PM UTC-5), the total size of that directory is 726.6236514812 GB that includes 1256 files in it. Thanks anyway.
Question: Hello,

I have a list of Uniprot/swissprot identifiers, such as:

P12281
Q05397
P50430
P54904
...

I also have the short names for the corresponding identifiers :

MOEA_ECOLI
FAK1_HUMAN
ARSB_RAT
P5CR1_ARATH
...

Does anyone knows how can I get the full protein names for these proteins in batch?, like the first one would be:
Molybdopterin molybdenumtransferase

If I could get even more information in this entry, like the associated GO terms, even better...

thank you!

Answer: Try the [Retrieve/ID mapping][1] tool from UniProt. MOEA_ECOLI FAK1_HUMAN ARSB_RAT P5CR1_ARATH will be your 'UniProtKB AC/ID'
and if you choose UniProtKB you will get a table with the protein names. You can filter the [results][2] by Gene Ontology, Taxonomy and others.

[1]: http://www.uniprot.org/mapping/
[2]: http://www.uniprot.org/uniprot/?query=yourlist:M2016052483C3DD8CE55183C76102DC5D3A26728B9DFB5AV&sort=yourlist:M2016052483C3DD8CE55183C76102DC5D3A26728B9DFB5AV&columns=yourlist(M2016052483C3DD8CE55183C76102DC5D3A26728B9DFB5AV),id,entry%20name,reviewed,protein%20names,genes,organism,length
Question: hi,

i have a correlation matrix of genes but i am going to remove alphabet from rows and columns to prepare my file for matlb i did like below

> mycounts <- read.table("tmatnorm_rld500.txt", header = T, sep = "\t", row.names=1)
> head(mycounts[,1:3])
AT1G01060 AT1G01170 AT1G01180
AT1G01060 1.0000000 0.3193360 0.6272994
AT1G01170 0.3193360 1.0000000 0.3178880
AT1G01180 0.6272994 0.3178880 1.0000000
AT1G01260 0.2658624 0.3588030 0.2542234
AT1G01380 0.6178751 0.1561297 0.7052692
AT1G01490 0.5990499 0.6129301 0.6424225
> rownames(mycounts) <- NULL
> head(mycounts[,1:3])
AT1G01060 AT1G01170 AT1G01180
[1,] 1.0000000 0.3193360 0.6272994
[2,] 0.3193360 1.0000000 0.3178880
[3,] 0.6272994 0.3178880 1.0000000
[4,] 0.2658624 0.3588030 0.2542234
[5,] 0.6178751 0.1561297 0.7052692
[6,] 0.5990499 0.6129301 0.6424225
> colnames(mycounts) <- NULL
>
> head(mycounts[,1:3])
[,1] [,2] [,3]
[1,] 1.0000000 0.3193360 0.6272994
[2,] 0.3193360 1.0000000 0.3178880
[3,] 0.6272994 0.3178880 1.0000000
[4,] 0.2658624 0.3588030 0.2542234
[5,] 0.6178751 0.1561297 0.7052692
[6,] 0.5990499 0.6129301 0.6424225

> write.table(mycounts, file = "correlation.txt", dec = ".", sep = "\t", quote = FALSE, row.names=F)
> mycounts <- read.table("correlation.txt", header = T, sep = "\t")
> head(mycounts[,1:3])
V1 V2 V3 <---this line
1 1.0000000 0.3193360 0.6272994
2 0.3193360 1.0000000 0.3178880
3 0.6272994 0.3178880 1.0000000
4 0.2658624 0.3588030 0.2542234
5 0.6178751 0.1561297 0.7052692
6 0.5990499 0.6129301 0.6424225

but as you consider, i have something in colnames?? then how I can get rid of alphabet and only have the similarity values please?

thank you

Answer: write.table(m,file="outfile,txt",sep="\t", col.names = F, row.names = F)
Question: I have been asked to recommend introductory books and resources to R and Bioconductor. My problem is just, I never read a book to learn R or Bioconductor, so I have no experience with this and cannot recommend one. I am interested in mainly introductory books, possibly targeting various groups of readers (computer scientists, molecular biologists, (bio-)statisticians), any recommendation appreciated.

For example, I used the following resources:

The R-manuals, especially the R intro

There are also a lot of contributed documents there on the R web site, but I didn't
use.

If a package from Bioconductor interests me, I read the package vignette.

I read the Bioconductor mailing list, that helps to see what other people use.

I have the "Venebles, Ripley. S Programming" book, that is hardly introductory.

Which books did you find helpful or completely useless to learn R/Bioconductor? For example: R Programming for Bioinformatics looks promising, anybody read it?

Or do you share my reluctance towards R-books and prefer online resources?

Answer: Oh, I'm glad this got bumped--just yesterday I saw a new guide to R that might be useful:

The Undergraduate Guide to R
Question: I am trying to filter multi-allelic sites in a VCF file to find the true multi allelic sites. Some of these have 1 read in 200, that I want to get rid of. The other alternative looks good. However, when I try to filter the alleles lower than frequency (0.05) with `vcftools` it gets rid of all variant not just the minor allele. Is there a way to filter out the minor alternate according to the frequency but not the variant completely?

Answer: Found the answer!

1. Convert multiallelic to biallelic vcf first

bcftools norm -m - file.vcf > biallellic.vcf

2. Filter the alternative alleles under certain value

bcftools view -e "FORMAT/AD[:1]<2 && INFO/AD[1]<5" biallelic.vcf > biallelic-filtered.vcf

3. Convert biallelic vcf to multiallelic vcf

bcftools norm -m + biallelic-filtered.vcf > multiallellic-filtered.vcf
Question: I currently have a heatmap (which is k-means clustered) of the top 1000 most variable genes, but I want to be able to narrow it down to about 100-200. The reason why it's set so high is that when I want to analyse each cluster, I end up filtering out any genes that are not considered significantly differentially expressed genes (e.g. padj < 0.05, log2foldchange < 0 or > 0). However when I do this, I get left with a small amount of genes which isn't ideal for when I'm doing GO-term analysis - hence why I start off with a larger number (like 1000), so that even after filtering there's enough genes in the cluster for me to analyse.

I was wondering if there is a way to filter out the genes beforehand, so that all the genes in the heatmap are already significantly expressed and I won't have to filter them out afterwards - that is, all the genes in the heatmap will actually be used for downstream analyses. Hope that makes sense! Also open to any recommendations and suggestions on how to improve my approach.

Here is my code:

#heatmap
topVarGenes <- head(order(rowVars(assay(rld)), decreasing = TRUE), 1000)

#clustered heatmap
set.seed(1234)

k <- pheatmap(assay(rld)[topVarGenes,], scale="row",kmeans_k = 4)
clusterDF <- as.data.frame(factor(k$kmeans$cluster))
colnames(clusterDF) <- "Cluster"

OrderByCluster <- assay(rld)[topVarGenes,][order(clusterDF$Cluster),]

pheatmap(OrderByCluster,
scale="row",annotation_row = clusterDF,
show_rownames = FALSE,cluster_rows = FALSE)
acute_hi <- rownames(clusterDF[clusterDF$Cluster == 1,,drop=FALSE])

#DEGs
resSigind = res[ which(res$padj < 0.05 & res$log2FoldChange > 0), ]
resSigrep = res[ which(res$padj < 0.05 & res$log2FoldChange < 0), ]
resSig = rbind(resSigind, resSigrep)

#filtering out any genes from the heatmap which are not significantly different (e.g. do not overlap with DEGs)
acutehi_cluster <- resSig[acute_hi,]

Answer: Hi,

What you can do is to retrieve the list of genes that are differentially expressed, and use this list first to filter the `rlog` transformed matrix that you have, before trying to obtain the most variable genes. This would return the most variable genes among the genes that are differentially expressed.

I hope this answers your question,

António
Question: Hi,

I was asked to create a dashboard hosting some genome-wide microarray and RNA-seq results that were computed using limma (an empirical Bayes method). Basically, the dashboard would allow users to retrieve gene(s) and see if it was differentially expressed (DE) or not, mainly to generate new hypotheses.

In the past, when I wanted to test a specific subset of genes (e.g, pathway X, ~50 genes) from a genome-wide experiment, I'd proceed as usual, but in the last step, I'd only adjust for multiple testing the subset list of genes of interest and not even look at the remaning genes (as they were not of interest).

Thinking in this dashboard, please consider this scenario:

User inputs 100 genes in the dashboard to screen (representing 1 or 2 solid biological pathways). Then, the original DE table is filtered and Benjamini-Hochberg or Benjamini-Yekutieli FDR are applied only to those 100 genes and returned to user. User see the results and then decides to include more genes from the same or other pathways, now totalling 150 genes. **What should the dashboard do?**

- A) Adjust P-values for all 150 genes, not the new 50;
- B) Adjust P-values separately for each set (adjust for 100) and then adjust for 50;
- C) Best to not use genome-wide methods (e.g, empirical Bayes), perform DE with standard NHST (e.g, t-test, ANOVA) and adjust for multiple testing using 150 genes;
- D) Best to not use genome-wide methods (e.g, empirical Bayes), perform DE withstandard NHST (e.g, t-test, ANOVA) and adjust for multiple testing separately for set 1 and 2;
- E) Other better appraches?

I also know that these procedures do not adjust confidence intervals, so that's another thing...

Does it make sense? I have a feeling that this is not ideal. But I'm not sure if there's a way to proper control the error rates using this dashboard. I also assume the user would be honest and tell the software which lists were already retrieved or the system would have some sort of caching and editing past queries.

Probably my question does not has an easy answer, but I'd like to hear from statisticians and more experienced peers.

Thank you.

Answer: It's generally not advised to run FDR correction on subsets of RNA-seq values. FDR uses the distribution of p-values for correction; importantly the somewhat uniform right tail of the distribution is used to find the cutoff value. Depending on how you are sampling p-values, you are changing the distribution of p-values going into the calculation, and these distributions may not be representative of the whole data.

A better approach would be to apply FDR correction on the whole dataset, and simply subset that precalculated table on your dashboard.
Question: Hi Everyone,

I am new to phylogenetics and using RAXML to generate my phylogenetic trees. I did bootstrapping (1000) to get the support values on the best ML tree. I get two output files: one is the ML best tree and the other one contains the bootstrap values. Which tree I should use as an output? They both look pretty similar.

Thanks, RT

Answer: They will be identical, except the bootstrapped tree will have support values (generated through bootstrap resampling).

For systematic inference, a tree without support values is more-or-less meaningless, as you might have bipartitions that are not well-supported yet still present in your consensus tree.

For comparative phylogenetic methods, you ought to 1) use the best-supported tree or 2) use a distribution of trees.
Question: I have a de novo transcriptome for an ant species of about 100k transcripts assembled by Trinity. I've quantified expression using [Sailfish](http://www.cs.cmu.edu/~ckingsf/software/sailfish/) and performed some analyses to identify genes with significant expression patterns.

Then, I became aware of the program [DeconSeq](http://deconseq.sourceforge.net/) to remove contaminants. I ran the program and removed ~5k of my transcripts, and spot-checking showed they all BLAST to bacteria or human. All good so far. Then, I re-ran Sailfish on the "cleaned" transcriptome of about 99k transcripts. My naive expectation was that this should barely effect the results - the reads that mapped to the "contaminants" shouldn't map at all. Instead, I find substantial changes, with about twice as many genes having significant expression patterns. I checked the unmapped ratio and it stays the same between mappings to the two transcriptome files.

So - my conundrum is how to deal with mapping reads to known contaminants?

Do I:

a) Map reads to the complete transcriptome, including contaminants, and then remove the known contaminants. In this case, it seems that I risk incorrectly mapping reads to the contaminants and losing information (false negatives), analogous to a Type I (producer's) error.

b) Map reads to the "cleaned" transcriptome. This seems analogous to a Type II (consumer's) error where I risk finding significant changes in expression that are due to incorrectly mapping 'contaminant' reads to true transcripts (false positives).

Any thoughts appreciated!

cross-posted on Sailfish user group [here](https://groups.google.com/forum/?utm_medium=email&utm_source=footer#!msg/sailfish-users/bCzPIg3fTps/8gYOfjZlhiEJ)

Answer: 1. "Mapping reads to contaminant sequences" - If the reads are without contamination, obviously the contaminant sequences will not get mapped. You can safely remove the contaminated assembled sequences by putting a coverage cutoff.
2. "Incorrectly mapping 'contaminant' reads to true transcripts" - Again, if the reads are cleaned (without contamination), then is no need to worry about expression level of true transcripts.

Now, if I were you, even before mapping the reads to find the expression levels, I would remove the contaminant reads by mapping them against contaminant (Bacterial, Viral genome or rRNA) Sequences. To confirm they are contaminant reads, If possible, I would take reference genome or nearest reference species sequence and map these contaminant reads against them. If the reads map to both contamination and nearest species, I would retain to map to my transcriptome, but I would throw away those which hit only Contamination. Even though, one may remove contaminants from assembly, if assembled with contaminant reads, there is always a chance of misassembly.
Question: Hello there, I have 2 files with the same header and same 1st column, I would like to compare rows by rows and calculate
the percentage of matches of the genotype. Somehow I try to transpose it to compare column by column by using the function compare in r, I even try use awk to do the comparing, somehow I couldn't make it. Anyone please help me with some hints or tips. Thanks

**Set1**

barcode	SNP000072	SNP000119	SNP000179	SNP001106	SNP001150
165974-1	A:A	A:A	A:A	G:G	C:C
165974-2	A:A	A:A	A:A	G:A	C:C
165974-3	A:A	A:A	G:A	G:A	C:C
165974-4	A:A	A:A	A:A	A:A	C:A
165974-5	A:A	A:C	A:A	G:A	?

**Set2**

barcode	SNP000072	SNP000119	SNP000179	SNP001106	SNP001150
165974-1	A:A	A:A	A:A	G:G	C:C
165974-2	A:A	A:A	A:A	G:A	C:C
165974-3	A:A	A:A	A:A	G:A	C:C
165974-4	A:A	A:A	A:A	A:A	C:A
165974-5	A:A	A:A	A:A	G:A	C:C

**Expected output**

barcode percentage(%)
165974-1	100
165974-2	100
165974-3	80
165974-4 100
165974-5 60

Answer: join -t $'\t' -1 1 -2 1 \
<(sort -t $'\t' -k1,1 input1.txt)\
<(sort -t $'\t' -k1,1 input2.txt) |\
grep -v barcode |\
awk -F '\t' 'BEGIN{printf("barcode\tpercent\n");}{NSAMPLES=(NF-1)/2;T=0.0;for(i=2;i<2+NSAMPLES;i++) {j=i+NSAMPLES; if($i==$j) T++;} printf("%s\t%d\n",$1,100*(T/NSAMPLES));}'

barcode	percent
165974-1	100
165974-2	100
165974-3	80
165974-4	100
165974-5	60
Question: Hello all

I have some gene IDs from Ensembl and I want to get their transcripts' exons and introns (sequences) so that later I can determine exon/intron boundaries and do some analyses.

I've discovered Ensembl Rest API, which is a really easy and clean way of getting data and played around with it a bit. Using this API, I could get coding transcripts of the genes and then sequences of these transcripts. However, I couldn't find any way to distinguish exonic and intronic regions in these sequences.

Here is my script that gets sequences of transcripts of "ENSG00000197568" gene in FASTA format. And I want to get exons and introns like Ensembl gives us in [here][1].

[1]: http://www.ensembl.org/Homo_sapiens/Transcript/Exons?db=core;g=ENSG00000197568;r=1:70820488-70851022;t=ENST00000359875

```
#!/usr/local/bin/python

import httplib2, sys, re, json

def check_response(response):
if not response.status == 200:
print "Invalid response: ", response.status
sys.exit()

http = httplib2.Http(".cache")
server = "http://beta.rest.ensembl.org/"
gene_id = "ENSG00000197568"
query = "sequence/id/" + gene_id + "?type=cds;multiple_sequences=1"
content_type = "text/x-fasta"

response, content = http.request(server + query, method="GET", headers={"Content-Type":content_type})
check_response(response)

transcripts = re.findall(">(.*)", content)
f = open("output.fasta", "wa")

for transcript in transcripts:
query = "sequence/id/" + transcript
response, content = http.request(server + query, method="GET", headers={"Content-Type":content_type})
check_response(response)
f.write(content)

f.close()
```

Thanks in advance

Answer: Hi Gungor,

We have taken your suggestion into account and added an option to softmask intronic regions

This option is available on our new rest server, http://rest.ensembl.org, along with improved performance.

The following: http://rest.ensembl.org/sequence/id/ENSG00000157764?content-type=text/plain;mask_feature=1

will return the whole gene sequence, with intron sequences in lower case.

In your example, http://rest.ensembl.org/sequence/id/ENSG00000197568?content-type=text/plain;type=cds;multiple_sequences=1

you are already retrieving only the coding sequence for each transcript in the gene. Hence, there are no intronic regions.

I hope this helps and please do not hesitate to contact us if you have any further enquiries.

Regards,
Magali
Question: I want to convert a paired SRA file to FASTA file to do some downstream analysis. For example the following link that is a metagenomics wgs paired sequencing:
http://www.ncbi.nlm.nih.gov/sra/?term=ERR1018199
I know I can use fastq-dump and I did same as bellow:

./fastq-dump myfile.sra --fasta

and it gives me:

>ERR1018199.1 FCD0JMKACXX:8:1101:1051:3172#ACTTGAAT length=200
GGAAGCGGTGTTCTGTTTCTCCTTCCATATATTTCCACAGACTGGCATATTCTTCTTTCTTCTCGCCCAT
TTCCTGCAATTGCAGCATATATTGCCAGATGATGAAAAAAATATACATCAAGGAAAACAAGTCTATCTGG
CAATATATGCTGCAATTGCAGGAAATGGGCGAGAAGAAAGAAGAATATGCCAGTCTGTGG
>ERR1018199.2 FCD0JMKACXX:8:1101:1348:3168#ACTTGAAT length=200
CTTGCAGATTCTACAAAAAGAGTGTTTCATAAACTGGTCTATCAAAAGAAAGGTTAAACTCAGTGAGTTG
AACCCACACATCACAAAGTAGCTTCTGAGATATGTGGGTAATATCTGTATGGATGTTTGTATGATTGATG
TTACTGACATTGATTGCAAAGAAGGCGACAGCGTTGAGATTTTCGGAGATCATCTTCCTA
>ERR1018199.3 FCD0JMKACXX:8:1101:1451:3174#ACTTGAAT length=199
GTATCATGACCGGTCGTTCGGGCAACAACATTTGGTGTATCAGTCCGATGTTCGACCTCAACAAACCGAC

and if I use --split-files option. It gives me:

ERR1018199_1.fasta:

>ERR1018199.1 FCD0JMKACXX:8:1101:1051:3172#ACTTGAAT length=100
GGAAGCGGTGTTCTGTTTCTCCTTCCATATATTTCCACAGACTGGCATATTCTTCTTTCTTCTCGCCCAT
TTCCTGCAATTGCAGCATATATTGCCAGAT
>ERR1018199.2 FCD0JMKACXX:8:1101:1348:3168#ACTTGAAT length=100
CTTGCAGATTCTACAAAAAGAGTGTTTCATAAACTGGTCTATCAAAAGAAAGGTTAAACTCAGTGAGTTG
AACCCACACATCACAAAGTAGCTTCTGAGA
>ERR1018199.3 FCD0JMKACXX:8:1101:1451:3174#ACTTGAAT length=100
GTATCATGACCGGTCGTTCGGGCAACAACATTTGGTGTATCAGTCCGATGTTCGACCTCAACAAACCGAC
TATTTTACTAAATTCCCACATCGATACCGT
>ERR1018199.4 FCD0JMKACXX:8:1101:1254:3223#ACTTGAAT length=100

ERR1018199_2.fasta:

>ERR1018199.1 FCD0JMKACXX:8:1101:1051:3172#ACTTGAAT length=100
GATGAAAAAAATATACATCAAGGAAAACAAGTCTATCTGGCAATATATGCTGCAATTGCAGGAAATGGGC
GAGAAGAAAGAAGAATATGCCAGTCTGTGG
>ERR1018199.2 FCD0JMKACXX:8:1101:1348:3168#ACTTGAAT length=100
TATGTGGGTAATATCTGTATGGATGTTTGTATGATTGATGTTACTGACATTGATTGCAAAGAAGGCGACA
GCGTTGAGATTTTCGGAGATCATCTTCCTA
>ERR1018199.3 FCD0JMKACXX:8:1101:1451:3174#ACTTGAAT length=99
GGGAGGTAGTTGGGGCAATACGCTTTCTATTCCGTTTTTTCCGGAAACTTCTTCTTCACATGAGGCAAGA
AAGATCAGATTGTATGCCTGCTCGGTTAT
>ERR1018199.4 FCD0JMKACXX:8:1101:1254:3223#ACTTGAAT length=100

My questions is, the first command simply concatenate the two reads with same read number. is it true? Can I use that for further analysis (for example assemble) or I have to use --split-files to produce two different files and apply my pipeline on those files?

Or more simpler can we say that we have the following sequence in our (meta)genome? or not?

>ERR1018199.1 FCD0JMKACXX:8:1101:1051:3172#ACTTGAAT length=200
GGAAGCGGTGTTCTGTTTCTCCTTCCATATATTTCCACAGACTGGCATATTCTTCTTTCTTCTCGCCCAT
TTCCTGCAATTGCAGCATATATTGCCAGATGATGAAAAAAATATACATCAAGGAAAACAAGTCTATCTGG
CAATATATGCTGCAATTGCAGGAAATGGGCGAGAAGAAAGAAGAATATGCCAGTCTGTGG

Answer: The concatenated file is mainly for convenience, for example to move the data around as a single file. In general it cannot be used for downstream analysis.

There is meaning for the data in each pair and most tools will be unable to split the data unless it comes from different files.