Use of Sine Shaped High-Frequency Rhythmic Visual Stimuli Patterns for SSVEP Response Analysis and Fatigue Rate Evaluation in Normal Subjects

Background: Recent EEG-SSVEP signal based BCI studies have used high frequency square pulse visual stimuli to reduce subjective fatigue. However, the effect of total harmonic distortion (THD) has not been considered. Compared to CRT and LCD monitors, LED screen displays high-frequency wave with better refresh rate. In this study, we present high frequency sine wave simple and rhythmic patterns with low THD rate by LED to analyze SSVEP responses and evaluate subjective fatigue in normal subjects. Materials and Methods: We used patterns of 3-sequence high-frequency sine waves (25, 30, and 35 Hz) to design our visual stimuli. Nine stimuli patterns, 3 simple (repetition of each of above 3 frequencies e.g., P25-25-25) and 6 rhythmic (all of the frequencies in 6 different sequences e.g., P25-30-35) were chosen. A hardware setup with low THD rate (<0.1%) was designed to present these patterns on LED. Twenty two normal subjects (aged 23–30 (25 ± 2.1) yrs) were enrolled. Visual analog scale (VAS) was used for subjective fatigue evaluation after presentation of each stimulus pattern. PSD, CCA, and LASSO methods were employed to analyze SSVEP responses. The data including SSVEP features and fatigue rate for different visual stimuli patterns were statistically evaluated. Results: All 9 visual stimuli patterns elicited SSVEP responses. Overall, obtained accuracy rates were 88.35% for PSD and > 90% for CCA and LASSO (for TWs > 1 s). High frequency rhythmic patterns group with low THD rate showed higher accuracy rate (99.24%) than simple patterns group (98.48%). Repeated measure ANOVA showed significant difference between rhythmic pattern features (P < 0.0005). Overall, there was no significant difference between the VAS of rhythmic [3.85 ± 2.13] compared to the simple patterns group [3.96 ± 2.21], (P = 0.63). Rhythmic group had lower within group VAS variation (min = P25-30-35 [2.90 ± 2.45], max = P35-25-30 [4.81 ± 2.65]) as well as least individual pattern VAS (P25-30-35). Discussion and Conclusion: Overall, rhythmic and simple pattern groups had higher and similar accuracy rates. Rhythmic stimuli patterns showed insignificantly lower fatigue rate than simple patterns. We conclude that both rhythmic and simple visual high frequency sine wave stimuli require further research for human subject SSVEP-BCI studies.

EEG-steady state visual evoked potential (SSVEP) signal is one of the most promising modalities that has several comparative advantages including high signal to noise ratio (SNR) (Ahn et al., 2016), high information transfer rate (ITR) and minimal requirement for training of subject (Middendorf et al., 2000;Wu et al., 2008;Guger et al., 2012). SSVEP signal is a natural response of brain to periodic visual stimulation in range of 1-90 Hz (Herrmann, 2001). It usually includes a sinusoidal waveform with the same fundamental frequency as of the external visual stimulus, harmonics and occasionally subharmonics, and is mostly observed in occipital and parietal lobes (Zhang et al., 2012a).
In the SSVEP-based BCI, an extensive variety of frequencies have been utilized for visual stimulus generation and it is observed that SSVEP has the strongest amplitude when the flickering frequency is about 15 Hz (Pastor et al., 2003). Stimulus frequencies in the medium (12-25 Hz) and low frequency (up to 12 Hz) ranges give better SNR in SSVEP. However, these are associated with subjective fatigue and discomfort and higher risk of photosensitive epileptic seizures (Pastor et al., 2003;Duszyk et al., 2014). On the other hand, high-frequency (25-50 Hz) stimulus range yield lower amplitude in SSVEP (Diez et al., 2011;Ajami et al., 2018); however, they give stable subjective performance with the passage of time making them suitable for SSVEP-based BCI (Won et al., 2015). Several studies have utilized low to medium range frequencies as these have higher SNR than high frequency range and are easy to detect (Friman et al., 2007;Muller et al., 2008;Ortner et al., 2011). Recent studies have detected SSVEP by using visual stimuli at the range of 30-50 Hz (Zhu et al., 2010;F. Zhang et al., 2015). For example, controlling robotic wheelchair and computer mouse at 37-40 Hz frequencies (Diez et al., 2011(Diez et al., , 2013 and using 30-39 Hz frequencies for an optimized simple 7-letter speller (Chabuda et al., 2018).
Several studies have developed methods for the detection of SSVEP at high frequencies. For example, spatial filtering approach has been used for the SSVEP detection at 30, 35, 40, and 45 Hz (Molina and Mihajlovic, 2010). In addition, empirical mode decomposition (EMD) method has also been used for the SSVEP recognition at 25, 30, 35 Hz (Liu et al., 2014;Zhang et al., 2015).
The information transfer rate of a BCI system is proportional to the number of commands that are simultaneously available for the user (Yuan et al., 2013). In SSVEP paradigm, a narrow band of frequencies is available so as to have low eye fatigue and acceptable SNR. On the other hand, the refreshing rate of monitors is an important limiting factor for the number of available visual stimuli. Therefore, the number of available targeted frequencies for the modulation is an important problem that remains to be solved for benchmark BCI application such as spellers (Zhang et al., 2012b). Some studies have attempted to overcome this problem by using phase information to code more targets Jia et al., 2011) or by embedding more LEDs with a low precise frequency resolution interval between the targeted frequencies (Hwang et al., 2012). Others have proposed the sum of two or more frequencies applied simultaneously to a visual stimulus, mostly with square waveform (Srihari Mukesh et al., 2005;Shyu et al., 2010;Hwang et al., 2013). However, as reported in several studies, the use of square waveform results in appearance of odd harmonics in the SSVEP response that can negatively affect the detection accuracy and increases the ambiguity in choosing the best stimulus frequencies for the target recognition (Srihari Mukesh et al., 2005;Shyu et al., 2010;Zhang et al., 2012bZhang et al., , 2015Zhao et al., 2017).
In addition, the above mentioned limitations have been addressed by considering specific aspects of challenges. For example, use of high frequency for the reduction of fatigue (Sakurada et al., 2015); however, it decreases the BCI performance (Volosyak et al., 2011). As it is important to reduce fatigue rate and also have high accuracy rate for a real BCI application, several studies have used high frequency for fatigue reduction and low frequencies for improving the accuracy rate by introducing amplitude and frequency modulated visual stimulation methods (Chang et al., 2014;Dreyer and Herrmann, 2015;Dreyer et al., 2017). In these studies, a high frequency carrier was modulated by another frequency and the SSVEP was evoked at the difference of these two frequencies that covered low frequency range for improving the accuracy rate. This method reduced subjective perceptibility as well as fatigue rate; however, the accuracy rate in these studies was lower than the simple stimulus used in recent studies (Dreyer and Herrmann, 2015;Dreyer et al., 2017).
In this study, we used sequence of 3 different frequencies as high frequency rhythmic patterns for the reduction of subjective fatigue rate. We hypothesized that the accuracy of rhythmic patterns would be comparable with the recently used simple visual stimuli. Thus, we designed high frequency simple and rhythmic sine wave patterns emitted by LED with high precision and low total harmonic distortion (THD) rate. We also investigated the accuracy by using three different frequency detection algorithms. Finally, the discrimination of patterns was evaluated by the analysis of SSVEP responses and the subjective fatigue for the simple and rhythmic visual stimuli by using VAS.

Study Participants
Twenty-two healthy right handed (11 males and 11 females students; 23-30 years old, average age 25 ± 2.1 years) subjects with normal or corrected to normal vision and without any history of neurologic and psychiatric disorders and head trauma were enrolled in this study. They were recruited by the word of mouth or announcement via social media. All of them signed informed consent form based on the Declaration of Helsinki. The ethical approval for the study was given by the Medical Research Ethics Committee of the Tehran University of Medical Sciences.

Stimuli Design Paradigm
In the human visual system, three parallel pathways including Magnocellular (MC), Parvocellular (PC) and Koniocellular (KC), process and transfer visual information to visual cortex (Kaplan, 2014). Each pathway processes different physical parameters of visual stimuli characterized by their specific temporal and spatial resolutions (Purpura et al., 1990;Duszyk et al., 2014;Labecki et al., 2016). Magnocellular pathway is sensitive to difference in contrast and motion and depth information. Receptive fields of MC pathway neurons are larger than other pathways and are sensitive to quick transients in retinal stimulation. In other words, it is sensitive to high-frequency visual stimuli (Cheong et al., 2013;Duszyk et al., 2014). Information about red and green colors is carried by the PC pathway and it is also sensitive to the shape of stimulus and exhibits sustained response to retinal stimulation. Koniocellular pathway carries blue and yellow color information and responds to the spectrum of stimuli (Kaplan, 2014). Based on the above characteristics of human visual system, we expected high SSVEP generation by MC pathway because of its larger receptive fields, its relevant processing of visual stimuli and motion perception. The most common colors used for BCI applications are green, blue, and red and white (Zhu et al., 2010). For choosing the best color to have strong SSVEP response, we selected red color based on PC pathway sensitivity and for obtaining sustained SSVEP response (Zhu et al., 2010;Duszyk et al., 2014).
Knowledge about the shape and size of stimuli presenting device is also crucial for designing a BCI-paradigm. Firstly, we selected LED because of the limitation of the refresh rate of LCD (Zhao et al., 2017). The size of square shaped LED was 4 × 4 cm 2 . Rectangular, square or checkerboard are the most common shapes used in BCI studies and choosing any one of them, have no significant effects on response (Duszyk et al., 2014). Additionally, as increasing the size of LED evokes stronger SSVEP response, we kept 4 × 4 cm 2 size which was acceptable for use in SSVEP-BCIs (Duszyk et al., 2014). According to our paradigm, each high-frequency visual pattern was coded by the permutation of frequencies that are available in the [25-50 Hz] range. By choosing N frequencies from the frequency set, we could get N N different patterns. These patterns contained N same frequencies, (N N -N-N!) partly different and N! completely different frequencies. General block diagram for designing visual stimulus is shown in Figure 1.
For applying the above visual stimuli in our SSVEP study, three high frequencies 25, 30, and 35 Hz were selected to generate various patterns. We generated 9 patterns, (each having 3 frequencies) from 27 possible ones. Three patterns comprising of repetitions of same frequencies (P25-25-25, P30-30-30, and P35-35-35) served as the simple ones (Garcia, 2008). The other six patterns that served as rhythmic patterns, comprised of all of 3 frequencies in different sequences (P25-30-35, P25-35-30, P30-35-25, P30-25-35, P35-30-25, and P35-25-30). Use of the 3 frequencies in the 6 above sequences gave maximum frequency changes and unpredictability in rhythmic patterns. These 6 rhythmic patterns also have ascending (e.g., P25-30-35) descending (e.g., P35-30-25) and zigzag (e.g., P30-25-35) trends in frequencies. The order of presentation of simple and rhythmic visual stimuli pattern was random with 3 simple patterns interspersed between 6 rhythmic patterns and same for all the subjects as shown in Supplementary Table S1. Thus, we were able to investigate the effect of maximum frequency changes on subjective fatigue rate between simple and rhythmic group patterns. Break time of 2 s was placed between the trials for better adaptability. The time series diagram of visual stimuli presentation session for each simple and rhythmic patterns is shown in Figures 2, 3 respectively.

Stimulator Design
All patterns were first designed and entered in MATLAB software (Release 2016b, The MathWorks, Inc, Massachusetts, United States) with the sampling rate of 44100 sps and saved as .wav format to be sent to a custom-made DAC board and LED driver (Figures 4, 5). A precise DC bias was added to the sine waveform pattern with appropriate amplitude to lay in the linear region of operation of the LED. For the LED that was utilized in our setup, the linear region range was 0 to 60 milliamperes with 5.5 Volt DC bias. The stimulator hardware FIGURE 1 | Block diagram of designated high-frequency visual stimuli patterns used in the study. Colors in each segment indicate different frequencies. i,j,k indices relate to individual high frequencies that were used within a segment. Ti, Tj, and Tj show duration of representation of each fi, fj, and fk frequency respectively, in a segment and Trest expresses the rest time before presentaion of each single session. N is the total number of segments in the designed visual stimuli presentation session for a single pattern and shows the repetition number of a single segment. Tbreak denotes the break time between two consecutive segments. A trial contains the sequence of three frequencies fi,fj, and fk presentation. The duration of a visual stimuli presentation session for a single pattern was calculated by following equation: T session = T rest + (T Trial + T break ) × N.
FIGURE 2 | Six rhythmic visual stimuli presentation session. In each session, after a 10 s rest time, the selected pattern was presented as a trial for 6 s. Between the two consecutive trials, a 2 s break time was given. Therefore, each session that contained 10 trials took 90 s. At the end of each session, fatigue rate was reported by the subjects. The session finished after 2 min rest period.
was designed by the authors (MF., AM., BM., and AJ.) so that all the parameters were precisely adjusted, and finally the THD rate of ≤0.1% was achieved (see Supplementary Figures S1, S2). The National Instruments DAQ was used for recording optic sensor signal and the stimuli that were represented concurrently. A linear optical sensor (Texas Instruments) was used for recording the represented stimuli signal from LED (Figures 4, 5).
As shown in the Supplementary Figures S1, S2, the linearity of luminance was analyzed by placing the optical sensor connected to the stimulator hardware (specifications not mentioned due to the pending patent) in front of the LED. We produced a single  pure 20 Hz sinusoidal signal with the designed hardware. This signal was applied to the LED with a precision driver. The LED output (e.g., luminosity) was measured with the linear optical sensor. The sensor output was analyzed to calculate the THD of the stimuli.

EEG Data Acquisition and Subjective Fatigue Evaluation
Before starting the experiment, the data acquisition paradigm was explained to each subject. The experiment was carried out in a dimly lit room without electromagnetic shield, and the subjects sat on the comfortable chair at a distance of about 75 cm from the LED screen as shown in the Figures 4, 5. The order of visual stimuli presentation was same for all the subjects (see Supplementary Table S1). They were given 2 minrest after each session. Subjective fatigue rate was evaluated during the rest time by using a 10 point (0 = no fatigue and 10 = maximum fatigue) Visual Analogue Scale (VAS) (Shahid et al., 2011). The EEG signals were recorded continuously by a g.USB amplifier (g.tec, Graz, Austria) system using g.LADYbird active electrodes at Oz, O1, and O2 connected to 3 channels (Ch1, Ch2 and Ch3) for the reported strong SSVEP response. Right ear lobe and Fpz electrodes were selected as reference and ground respectively, by following the EEG international 10-20 standard system and all the 3 channels were sampled at 1200 Hz sampling rate (Figures 4, 5).

Pre-processing
Based on the specifications of the optic sensor and synchronizer gtec.'s pulse signals, the initial time point of EEG recording was marked by gtec.'s synchronized pulse and the onset time of the stimuli presentation was determined by the recorded optic sensor signal. Finally, the beginning time point of EEG signal with respect to the visual stimuli representation was calculated by the subtraction of optic sensor initial point from gtec.'s pulse onset. Therefore, each segment of EEG-SSVEP was separated. We checked the signals visually and removed trend of data with detrend command in the MATLAB software. Single trials and mean of 10 trials were calculated for each pattern in order to utilize them for further respective single and mean trials data analysis.

PSD Based Analysis
SSVEP signals were filtered by offline 6-order Butterworth bandpass filter with cut-off frequencies of 10-40 Hz, (Diez et al., 2011). The power spectral density was calculated at each 2 s rectangular window, then normalized at each visual stimulus frequencies (25,30,35,and Hz) as: Where ch is the number of channels, M is the total number of channels, f i is the visual stimuli frequencies, PSDssvep is the power spectral density of EEG-SSVEP signal and PSDrest is the power spectral density of resting state EEG signal. Finally, max P (f i ) utilized for recognition of target stimuli as: Canonical Correlation Analysis (CCA) CCA basically works on two random sets of data that may have underlying correlation. The validity, and effectiveness of CCA has been shown in Bin et al. (2009b) study. In our study, first set of data was EEG-SSVEP signal and the second data set was visual stimulus signals. The assumption for utilizing CCA feature extraction method for SSVEP detection is that the source of SSVEP signal (x) is the output of a linear system with visual stimulus (Y) as the input signal. Y can be decomposed in Fourier series of its harmonics as: Where f is the fundamental frequency of visual stimulus, S is the sampling rate and T is the number of sample points. CCA algorithm finds a pair linear combination for x and y as, X = ax and Y = by to maximize the correlation between two canonical variables {X, Y}, based on the following optimization problem: The ρ i coefficients utilized as the canonical correlation coefficients for detection of targeted visual stimulus frequency as: For more details about the CCA we referred to Lin et al. (2007).

Least Absolute Shrinkage and Selection Operator Analysis (LASSO)
LASSO as a method for recognition of SSVEP was proposed by Zhang et al. (2012a). This analysis has the capability of model selection and shrinkage estimation methodology. Due to its sparse approximation constraints, this method provides low variance and high interpretable solution for linear regression. By considering a standard linear regression model that y represents in the SSVEP response observations as: Where the size of y is n×1 vector, X is a n×p matrix that demonstrates stimuli frequencies and their harmonics and ε denote noise vector with zero mean and constant variance. The LASSO solves the below optimization problem for estimating sparseβ vector as: Where || .|| 1 and || .|| 2 represent l 1 norm and l 2 norm and λ is a penalty parameter that helps to achieve sparse solution forβ. The contribution degree (CD) of each stimulus frequency and its harmonics in EEG response was calculated as: Where M is the total number of EEG channels and K is the number of considered harmonics andβ ch ch. h is theβ = [β 1 1.1 , β 1 1.2 ,β 1 1.3 ,β 1 1.4 ,β 1 1.5 ,β 1 1.6 ,. . . ,β 3 3.5 ,β 3 3.6 ] respect to each channel and harmonics. Finally, the maximum CD i is selected as targeted stimuli as:

Evaluation Methods
Classification accuracy was utilized for evaluating the discriminability of the patterns proposed in this study for SSVEP-BCI applications. For the mean of 10 trials classification, mean of 10 trials of SSVEP response of each pattern was calculated for each subject. Therefore, the total number of trials (for all nine patterns) was equal to 198 (= 9 patterns × 22 subjects). For single trial classification, the total number of trials was equal to 1980 (=10 trials for each pattern × 9 patterns × 22 subjects). Classification accuracy was formulated as: Acc = Number of valid trials Total number of trials × 100 (11) Where valid trials are the trials that were correctly classified. Classification steps carried out according to the block diagram shown in Figures 6 (see also Supplementary Figure S3) are given below: Step 1: Target frequency in window 1. The EEG data in window 1 was input of the PSD, CCA or LASSO algorithm to determine the dominated frequency (fi) in the current window.
Step 2: Target frequency in window 2. The EEG data in window 2 was input of the frequency recognition algorithm to determine the dominated frequency (fj) in the current window.
Step 3: Target frequency window3. Similar to that of Step 1 and 2, the dominated frequency (fk) was estimated based on the EEG in window 3.
Step 4: Target pattern recognition. After getting fi, fj and fk, we had the ordered sequence fi-fj -fk. We then compared it with the respective pattern. If it matched with the predefined pattern, the output was taken as 1 or else it was taken as 0.

Statistical Analysis
Statistical analysis was performed with SPSS (Version 16.0. Chicago, SPSS Inc. IBM Corp, released 2011) to evaluate the fatigue rate of 9 patterns used in this study. Appropriate nonparametric statistical tests were employed when fatigue rate data set failed to provide criteria for normal distribution. Therefore, first the Friedman test with a significance level α = 0.05 was employed as non-parametric test to overall demonstrate any significant fatigue rate changes in all visual stimuli patterns as multiple conditions, then the post hoc analysis was done to compare and specifically show each pair of patterns fatigue rate changes. Wilcoxon signed-rank test was applied for pairwise comparisons. Bonferroni correction was used because several pairwise comparisons were conducted in the experiment to get the critical value, so that the test had P < 0.0014 to be significant. Paired t-test was used for comparison of the simple and rhythmic pattern groups fatigue rate (α = 0.05). ANOVA with repeated measures was used to compare the amplitude of CCA and LASSO coefficients of SSVEP responses (α = 0.05) for single trial analysis. Unless mentioned otherwise, the results are expressed as mean ± standard deviation (S.D.). We compared three targeted frequencies coefficient's amplitude and then examined the patterns amplitude changes for simple and rhythmic groups. Post hoc analysis has done by using Bonferroni correction (set α < 0.017) in order to check the specific significant differences between each pairs. The CCA and LASSO accuracy results for individual window lengths were analyzed by Wilcoxon signed-rank test with α = 0.05. We statistically compared simple and rhythmic group accuracy results of mean of 10 trials with Wilcoxon signed-rank test (α = 0.05) and single trial accuracy results were compared with paired t-test (α = 0.05).

Results for Mean of 10 Trials Analysis
First, the PSD-based analysis was performed at an average of 10 trials in 2 s rectangular window and the recognition accuracies for 22 subjects were calculated as shown in Figure 7. For the Grand Average of SSVEP responses and PSD for all nine patterns, see Supplementary Figure S4. Then, the CCA and LASSO analysis were done for 0.5 s window lengths. By increasing the length of window by 0.25 s, the analysis (CCA, LASSO) was repeated for mean of 10 trials at 0.75, 1, 1.25, 1.5, 1.75, and 2 s window lengths as shown in Tables 1, 2 and Figure 8. There was no statistically significant difference between the mean accuracy results of CCA and LASSO.

Single Trial Analysis Results
All the three methods were also used for single trial analysis. The PSD results for single trial WL of 2 s are presented in Table 3.   CCA and LASSO analysis results are shown in the Tables 4, 5 for single trial in different WL durations (0.5 to 2 s). Comparison of the simple and the rhythmic groups visual stimuli was done for single trial and the mean of 10 trials in CCA and LASSO analysis and the respective results are shown in the Table 6. For single trials results, higher accuracy rate was seen for the rhythmic patterns ( Table 6 and Supplementary Table  S5). Paired t-test showed that the single trial accuracy analysis between the simple and the rhythmic groups had statistically significant differences for 2 s WL in both CCA [t (21) = −5.689, P < 0.001] and LASSO [t (21) = −5.086, P < 0.001] as shown in Table 6 and Supplementary Table S5.
Repeated measures ANOVA with a Greenhouse-Geisser correction showed that the amplitude of coefficients differed statistically significant between the frequency targets in rhythmic group (CCA method: [Target = 25 Hz, F (1.210, 1595.981) =  Figure S4 that showed the

Fatigue Rate Statistical Results
Friedman test showed that there was statistically significant difference between fatigue rate of 9 visual stimuli patterns χ 2 (8) = 49.275, P = 0.001. Post hoc analysis with Wilcoxon signed-rank tests was conducted with application of Bonferroni correction resulting in significance level of P < 0.0014. Descriptive statistics is shown in Table 7 and post hoc analysis results are shown in  Tables 7, 8 and Supplementary Figure S6). For the individual patterns, the maximum fatigue rate (4.95 ± 2.57) and minimum (2.90 ± 2.27) fatigue rate corresponded to P25-25-25 and P25-30-35 patterns respectively. Fatigue rate between the first (P25-30-35) and last (P30-25-35) pattern was not statistically significant (Table 8 and Supplementary Figure S6). However, there was a statistically insignificant trend toward higher fatigue rates at the later time points (Supplementary Figures S6 and S7).

DISCUSSION
In this study, we have used high-frequency low THD sine wave with simple and rhythmic patterns and presented them as visual stimuli. Additionally, we analyzed SSVEP responses with three well known methods (PSD, CCA and LASSO) in order to discriminate patterns. We also evaluated the effects of simple and rhythmic patterns on subjective fatigue rate in normal subjects. We found that the rhythmic patterns had significantly higher accuracy rates in single trial and overall insignificant higher rates in the mean of 10 trials, and resulted in comparatively lesser though insignificant subjective fatigue.
In comparison to the two apparently similar studies that didn't examine subjective fatigue rate and have utilized square pulse shape stimuli type without reporting THD rate (Zhang et al., 2012b(Zhang et al., , 2015, we used sine shaped high frequency range stimuli and also examined the simple and rhythmic group accuracies as well as their subjective fatigue rate (see Supplementary Table S4 for detailed comparison).
Most SSVEP studies have utilized low-frequency visual stimuli range that yields high SNR results and resultant high accuracy achievement, mostly higher than 80% (Zhang et al., 2012b;Sakurada et al., 2013). However, at lower frequency ranges (∼8 Hz), the risk of photosensitive epileptic seizures and high fatigue rate have to be considered. These two main challenges encouraged researchers to carry out SSVEP studies at high frequency range. At high frequency range, because of comparatively lower SNR, accuracy obtained is also lower (between 70 and 85%) than that of low-frequency studies (Garcia, 2008). We have attempted to improve the accuracy of high frequency stimuli by examining the rhythmic patterns and compared them with simple ones.

Validation of the Rhythmic Visual Stimuli Patterns
In our study, six permutations were designed with selection of three different high frequencies (25,30,and 35 Hz) in 3sequence including ascending, descending and zigzag rhythms (Figures 2, 3, and Supplementary Tables S2, S3). Three simple patterns had same frequency in the sequence.
Comparison of simple and rhythmic patterns showed higher accuracy rate for rhythmic patterns group (Table 6 and see  Supplementary Table S5). Recently a high-frequency SSVEP study has reported accuracy rates of 90-100% employing improved design paradigm and use of synchronized averaging of EEG trials (Zhang et al., 2015). We designed visual stimulus patterns at high-frequency range and obtained similar results (>90%) (Tables 1, 2). Our results for single trials analysis are in the range of 73-74% (Tables 4, 5) and are comparable with (Dreyer et al., 2017) that reported 33-87%. Our single trial accuracy rate was also comparable with a recent low frequency study (Zhang et al., 2012b). Additionally, our mean of 10 trials accuracy rates for rhythmic and simple groups are better than another high frequency SSVEP study (Zhang et al., 2015; see Supplementary Table S4). Use of appropriate sequence coding increases the unpredictability of visual patterns and leads to discriminative SSVEP response. Statistical analysis confirmed the discrimination of patterns by examination of the amplitude of SSVEP features. As shown in Supplementary  Tables S2, S3, simple, ascending, descending and zigzag patterns coefficient vectors had significant differences with each other (see Supplementary Figure S5). This is also seen in different patterns illustrated by CCA coefficients matrix (see Supplementary  Figures S3, S5).

PSD, CCA, and LASSO Results Comparison
Most SSVEP studies have used PSD based method for SSVEP recognition and so we first tried to analyze data with this method. For 2 s rectangular window mean accuracy rate was obtained to be equal to 88.35%. However, at high-frequency range, obviously, this method was not what we expected because of low amplitude in frequency spectrum of SSVEP response compared to the lowfrequency range (Garcia, 2008). CCA and LASSO methods could improve SSVEP-BCI speed compared to PSD. By utilizing these two methods, more robust discriminative features in shorter TW of EEG data can be extracted (Zhang et al., 2012a). For sequence coding to design the visual stimulus patterns, it was important to analyze the data in shorter TWs. As shown in the CCA and LASSO results (Figure 8), for TW = 0.75 s, we were able to obtain acceptable accuracy rate of higher than 90%. Our results of comparison between CCA and LASSO at different TWs were not significant. However, at low TWs length (≤1.25 s), LASSO did show better performance compared to CCA while at TW > 1.25 s, CCA showed better results than LASSO. These results are consistent with a previous study that reported better performance for LASSO at low TWs duration [0.5-2.5 s], while for TW > 2.5 s, CCA had same and even better results (Zhang et al., 2012a).
However, compared to our study, they used low frequency, square shaped visual patterns stimuli and had 9 participants (our study had 22 participants). Our CCA and LASSO results point to the importance of shape of visual stimulus. Based on the knowledge of authors, studies carried out so far that have compared LASSO and CCA analysis of SSVEP, have not considered stimulus harmonic distortion (Zhang et al., 2012a,b). These and other studies have used square pulse shape stimulus and therefore the odd harmonics of fundamental frequency existed in their visual stimulus frequency spectrum. Thus, the SSVEP response that was recognized in these studies did not correspond only to the fundamental frequency. In our study, by designing the sine shape visual patterns with low THD it was possible to have SSVEP response elicited by fundamental frequency almost devoid of distortion. Thus, our CCA and LASSO analysis of the SSVEP responses were largely

Fatigue Rate Evaluation
Using high frequency range stimuli can reduce the subjective fatigue rate compared to the low frequency (Diez et al., 2011). Studies that have examined the fatigue rate included simple high frequency (Diez et al., 2011(Diez et al., , 2013 while others have used the advantages of high frequency range as carrier signal (Chang et al., 2014;Dreyer and Herrmann, 2015;Dreyer et al., 2017). This was done to decrease visual perceptibility and consequently decrease the fatigue rate. In our study fatigue rate was evaluated for the high frequency rhythmic and simple patterns. We used VAS as psychometric response scale for fatigue rate evaluation. Our results show that there was overall no significant difference between subjective fatigue rates of simple and rhythmic stimuli pattern groups, though the VAS was lower for the rhythmic patterns group. While significant fatigue rate change was obtained within both simple and rhythmic visual stimuli groups, however, within group fatigue rate variations in the rhythmic group patterns were lower than the simple group. In addition, rhythmic group pattern P25-30-35 had the least VAS score of our study among all the individual patterns. Maximum and minimum fatigue rates were observed for P25-25-25 and P35-35-35 patterns respectively, in the simple group, indicating that higher frequencies cause reduced subjective fatigue. Maximum and minimum fatigue scores in the rhythmic group was for P35-25-30 and P25-30-35 patterns respectively, indicating the importance of the order of frequencies presented in a sequence (see Supplementary Figure S6).
In this study, we observed higher but insignificant trend in VAS values over time indicating higher subjective fatigue from beginning to the end of experiment (Supplementary Figures S6, S7). As subjective fatigue was expected to increase over time in our experiment, compared to few other recent SSVEP studies (Diez et al., 2011(Diez et al., , 2013Dreyer and Herrmann, 2015;Dreyer et al., 2017) we kept acceptable rest time (2 min) to minimize it. However, the aim of our study was not to evaluate fatigue at various time points for each pattern, rather we were more interested in knowing the fatigue rate between simple and rhythmic groups as well as its variation within these two groups. Thus, three simple group patterns were randomly interspersed within six rhythmic group patterns. This order was kept same for all the subjects so as to get authentic data. As a result, the patterns of rhythmic group were distributed in the first, middle and in the last time points of each experiment with simple group in between them. Overall, this minimized the influence of time in VAS analysis for simple and rhythmic groups. The duration of our experiment (32 min) might be seen as longer, as seen by higher VAS scores with time, however, evaluation of fatigue between consecutive patterns was not aimed in our study. We suggest future experiments to evaluate the fatigue rate with individual patterns.

Limitations of the Study and Suggestions for the Future Works
For designing visual stimulus patterns, several limitations have to be considered. Perhaps the most important point is the number of frequencies that are selected. Although, there is a trade-off between accuracy and speed of BCI applications which is related to the length of each trial, as when the number of frequencies in each sequence increases, the speed of the BCI decreases. However, still, a high accuracy rate can be obtained.
In this study, we allocated equal time for each frequency in  one trial. However, for future works, it is suggested that time of each frequency can be optimized and set differently for improving the speed. Additionally, for increasing the speed, it could be better to detect the lowest number of each frequency repetition for SSVEP generation and then setting the least time required for placing each frequency in the specific pattern. For example, higher frequencies may need lesser time compared to the lower frequencies, because in Ti duration time, the number of periods that higher frequencies are repeated is more than lower frequencies. Thus, SSVEP response to higher frequencies may be generated faster than low frequencies.
For N individual pattern analysis, we suggest that the order of patterns should be randomized and selected from existing N! permutations. However, this approach will obviously require more time and experiments for evaluation of the all permutations. In addition, it will increase the subjective fatigue rate over time and may lead to significant effect on the responses. Finally, within/between subject effects and responses should also be examined.
For the analysis of SSVEP, various methods of analysis are reported recently (Liu et al., 2014;Cao et al., 2015;Zhang et al., 2015). We only used three conventional methods (CCA, LASSO and PSD) for analysis of SSVEP signals. In future works, these as well as other recent methods can be used.
In a real BCI application, various patterns are applied simultaneously. Therefore, for designing an online synchronous or even asynchronous real BCI based on employed rhythmic patterns, appropriate modifications in hardware set up must be considered together with the selection of accurate frequency value.
VAS measurement tool was used for subjective fatigue rate evaluation in our study. It's important to quantify fatigue rate more accurately and not only based on spectral analysis (Cao et al., 2014). For future studies, using an objective method for fatigue evaluation will be more valuable for consistency between SSVEP-BCI applications. For fatigue rate reduction, our suggestion is to use different colors along with rhythmic patterns as this leads to significant effect on SSVEP responses (Duszyk et al., 2014).

CONCLUSION
High frequency sine shaped rhythmic and simple patterns had low THD, and higher and similar accuracy rates. Rhythmic stimuli patterns showed insignificantly lower fatigue rate than simple patterns. We conclude that both rhythmic and simple visual high frequency sine wave stimuli can provide a base for further human subject SSVEP based BCI studies.

AUTHOR CONTRIBUTIONS
AK did experiments, acquired and analyzed the data, wrote first draft of manuscript. ZS acquired and analyzed data, helped in writing manuscript. MF designed the hardware and setup for the experiment. ES acquired and analyzed data, helped in writing manuscript. MH experimental design and data interpretation, Co-supervisor. AM Project Co-director, Co-supervisor, data interpretation, designed the hardware and setup for the experiment. BM Co-supervisor, data interpretation, designed the hardware for the experiment. AJ Original idea for study, Project Director, Principal Investigator, and Supervisor.

FUNDING
This study was supported by Tehran University of Medical Sciences under the project Extracting high frequency rhythmic steady state visual evoked potentials for utilizing in SSVEP-BCI applications with grant number 94-04-30-30878.