Error: Unable to compute SFS without -s option

[jkonv...@knights.ucf.edu]

Hi All,

I'm trying to obtain maximum likelihood values from the SFS of four populations. I've tried it two ways. Once with all six pairwise combinations, each with the suffix "_jointMAF" and once with one combined SFS titled "4pops_MSFS.obs".

Here's the command line I'm using:
./fsc26 -t 4pops.tpl -n 50000 -e 4pops.est -M -L 40 -d -q

It comes up with the error: "Unable to compute SFS without -s option"

When I add "-s0" it does one iteration and gives me a bunch of ".arp" files.

I'm not sure what I'm doing wrong. When I did this with three populations it worked fine and gave me maximum likelihood values. With four populations it is much more complicated. Any help would be greatly appreciated.

Thanks,
Johnny

[Laurent Excoffier]

Hi,

it may be that your input file is badly formatted. Can you make it accessible?

best

laurent

[jkonv...@knights.ucf.edu]

Hi,

I have six pairwise observation files:

ModelB_4pops_jointMAFpop1_0.obs
1 observations
d1_0 d1_1 d1_2 d1_3 d1_4 d1_5 d1_6 d1_7 d1_8 d1_9 d1_10 d1_11 d1_12 d1_13 d1_14 d1_15 d1_16 d1_17 d1_18 d1_19 d1_20
507 17 43 18 26 11 37 18 32 33 27 0 0 0 0 0 0 0 0 0 0

ModelB_4pops_jointMAFpop2_0.obs
1 observations
d2_0 d2_1 d2_2 d2_3 d2_4 d2_5 d2_6 d2_7 d2_8 d2_9 d2_10 d2_11 d2_12 d2_13 d2_14 d2_15 d2_16 d2_17 d2_18 d2_19 d2_20 d2_21 d2_22 d2_23 d2_24 d2_25 d2_26 d2_27 d2_28 d2_29 d2_30 d2_31 d2_32 d2_33 d2_34 d2_35 d2_36 d2_37 d2_38 d2_39 d2_40 d2_41 d2_42 d2_43 d2_44 d2_45 d2_46 d2_47 d2_48 d2_49 d2_50 d2_51 d2_52 d2_53 d2_54 d2_55 d2_56 d2_57 d2_58 d2_59 d2_60 d2_61 d2_62 d2_63 d2_64 d2_65 d2_66 d2_67 d2_68 d2_69 d2_70 d2_71 d2_72 d2_73 d2_74 d2_75 d2_76 d2_77 d2_78 d2_79 d2_80 d2_81 d2_82 d2_83 d2_84 d2_85 d2_86 d2_87 d2_88 d2_89 d2_90 d2_91 d2_92 d2_93 d2_94 d2_95 d2_96 d2_97 d2_98 d2_99 d2_100 d2_101 d2_102 d2_103 d2_104 d2_105 d2_106 d2_107 d2_108 d2_109 d2_110 d2_111 d2_112
65 9 71 97 15 15 12 6 13 24 33 19 34 21 36 9 25 17 16 12 9 9 0 11 7 10 8 8 10 5 11 5 7 7 5 10 9 12 21 56 124 42 95 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

ModelB_4pops_jointMAFpop3_0.obs
1 observations
d3_0 d3_1 d3_2 d3_3 d3_4
364 75 69 0 0

ModelB_4pops_jointMAFpop2_1.obs
1 observations
d2_0 d2_1 d2_2 d2_3 d2_4 d2_5 d2_6 d2_7 d2_8 d2_9 d2_10 d2_11 d2_12 d2_13 d2_14 d2_15 d2_16 d2_17 d2_18 d2_19 d2_20 d2_21 d2_22 d2_23 d2_24 d2_25 d2_26 d2_27 d2_28 d2_29 d2_30 d2_31 d2_32 d2_33 d2_34 d2_35 d2_36 d2_37 d2_38 d2_39 d2_40 d2_41 d2_42 d2_43 d2_44 d2_45 d2_46 d2_47 d2_48 d2_49 d2_50 d2_51 d2_52 d2_53 d2_54 d2_55 d2_56 d2_57 d2_58 d2_59 d2_60 d2_61 d2_62 d2_63 d2_64 d2_65 d2_66 d2_67 d2_68 d2_69 d2_70 d2_71 d2_72 d2_73 d2_74 d2_75 d2_76 d2_77 d2_78 d2_79 d2_80 d2_81 d2_82 d2_83 d2_84 d2_85 d2_86 d2_87 d2_88 d2_89 d2_90 d2_91 d2_92 d2_93 d2_94 d2_95 d2_96 d2_97 d2_98 d2_99 d2_100 d2_101 d2_102 d2_103 d2_104 d2_105 d2_106 d2_107 d2_108 d2_109 d2_110 d2_111 d2_112
455 17 28 18 15 20 14 14 22 25 12 9 22 12 12 14 10 13 34 72 173 0 2 1 0 1 1 0 0 2 1 0 20 0 0 0 0 0 1 0 0 20 6 3 1 2 2 2 7 9 5 3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

ModelB_4pops_jointMAFpop3_1.obs
1 observations
d3_0 d3_1 d3_2 d3_3 d3_4
662 42 35 0 0

ModelB_4pops_jointMAFpop3_2.obs
1 observations
d3_0 d3_1 d3_2 d3_3 d3_4
574 7 25 0 0

Here's my template file:
ModelB_4pops.tpl
//Number of population samples (demes)
4 populations to simulate
//Population effective sizes (number of genes)
N_POP0
N_POP1
N_POP2
N_POP3
//Sample sizes
21
10
56
2
//Growth rates
0
0
0
0
//Number of migration matrices
4
//Migration matrix 0
0.000 0.000 0.000
0.000 0.000 0.000
0.000 0.000 0.000
//Migration matrix 1
0.000 0.000 0.000
0.000 0.000 0.000
0.000 0.000 0.000
//Migration matrix 2
0.000 0.000 0.000
0.000 0.000 0.000
0.000 0.000 0.000
//Migration matrix 3
0.000 0.000 0.000
0.000 0.000 0.000
0.000 0.000 0.000
//historical event: time, source, sink, migrants, new deme size, growth rate, migr mat index
3 historical event
TDIV23 2 3 1 1 RES23 0 1
TDIV1_ANC23 1 3 1 RES13 0 2
TDIV01 0 1 1 N_POP1 0 3

//Number of independent loci [chromosomes]
1 0
//Per chromosome: Number of linkage blocks
1
//per block: Datatype, numm loci, rec rate and mut rate + optional parameters
FREQ 1 0 MUTRATE

and my estimation file:
ModelB_4pops.est
// Priors and rules file
// *********************

[PARAMETERS]
//#isInt? #name #dist.#min #max
//all Ns are in number of haploid individuals
1 N_POP0 unif 10 1e6 output
1 N_POP1 unif 10 1e6 output
1 N_POP2 unif 10 1e6 output
1 N_POP3 unif 10 1e6 output
1 N_ANCAll unif 10 1e6 output
1 N_ANC23 unif 10 1e6 output
1 TDIV23 unif 10 1e6 output
1 TPLUS23 unif 10 1e6 hide
1 TPLUS1_ANC23 unif 10 1e6 hide

[RULES]

[COMPLEX PARAMETERS]

1 TDIV1_ANC23 = TDIV23+TPLUS23 output
1 TDIV01 = TDIV_ANC23+TPLUS1_ANC23 output

0 RES23 = N_ANC23/N_POP3 hide
0 RES13 = N_POP1/N_ANC23 hide

Thanks for any help.
Johnny

[Laurent Excoffier]

Hi,

first check your generated .simparam file and make sure that your input files are read correctly... it looks as if it was not the case

Besides MUTRATE is not defined in the est file

L

[angeles...@gmail.com]

Hi,

I'm having the same problem here, i.e. no estimates and lots of arp files. I've tried setting the divergence times either as rules or as complex parameters, but either way, I'm having the same message.

I'm running it as

./fsc26 -t NoMig.tpl -e NoMig.est -M -m --multiSFS -n10000 -N20000 -0 -l10 -L15 -C10 -c1 

This is the TPL file:

//Number of population samples (demes) 0-bairdi, 1-Phaeonotus, 2-caniceps,3-oreganus

4 populations to simulate
//Population effective sizes (number of genes)
N_POP0
N_POP1
N_POP2
N_POP3
//Sample sizes

16
16
16
16

//Growth rates
0
0
0

//Number of migration matrices
1

//Migration matrix 0
0.000 0.000 0.000 0.000
0.000 0.000 0.000 0.000

0.000 0.000 0.000 0.000
0.000 0.000 0.000 0.000
//historical event: time, source, sink, migrants, new deme size, growth rate, migr mat index
3 historical event

TDIV32 3 2 1 RELANC32 0 0
TDIV21 2 1 1 RELANC21 0 0
TDIV0_ANC21 0 1 1 RELANC01 0 0

//Number of independent loci [chromosomes]
1 0
//Per chromosome: Number of linkage blocks
1
//per block: Datatype, numm loci, rec rate and mut rate + optional parameters

FREQ 1 0 1.0e-8 OUTEXP

And here the EST file:

// Priors and rules file
// *********************

[PARAMETERS]
//#isInt? #name   #dist.#min  #max
//all Ns are in number of haploid individuals

1  N_POP0       unif     10    2e6   output
1  N_POP1       unif     10    2e6   output
1  N_POP2       unif     10    2e6   output
1  N_POP3       unif     10    2e6   output
1  N_ANCAll     unif     10    2e6   output
1  N_ANC21      unif     10    2e6   output
1  N_ANC32      unif     10    2e6   output
1  TDIV0_ANC21  unif     10    8e4   output
1  TDIV21       unif     10    8e4   output
1  TDIV32       unif     10    8e4   output
1  TPLUS       unif     10    8e4   output
1  TPLUS32       unif     10    8e4   output

[RULES]
//TDIV21 > TDIV32
//TDIV0_ANC21 > TDIV21

[COMPLEX PARAMETERS]

1 TDIV21 = TDIV32 + TPLUS32 output
1 TDIV0_ANC21 = TDIV21 + TPLUS output
0  RELANC32 = N_ANC32/N_POP2   hide
0  RELANC21 = N_ANC21/N_POP1    hide
0  RELANC01 = N_ANCAll/N_ANC21      hide

Many thanks in advance!

[Laurent Excoffier]

Same problem = same answer as above

Hi,

first check your generated .simparam file and make sure that your input files are read correctly... it looks as if it was not the case

L

You also 3 extra blank lines in your tpl file that are probably causing problems...

best 

L